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ABSTRACT:  The  CE-QUAL-ICM  three-dimensional  eutrophication  model  was  applied  to  Lake  Wash¬ 
ington,  northwestern  Washington  State,  for  the  period  1995-1997.  Transport  processes  were  obtained  from 
the  companion  CH3D-WES  hydrodynamic  model.  The  model  activated  18  state  variables  in  the  water 
column  including  physical  variables;  phytoplankton;  multiple  forms  of  carbon,  nitrogen,  and  phosphorus; 
dissolved  oxygen;  and  fecal  coliform.  The  model  provided  excellent  representation  of  the  annual  cycle  of 
temperature,  chlorophyll,  nutrients,  and  dissolved  oxygen  in  the  lake.  The  vertical  structure  was  also  well- 
represented.  The  water  column  was  coupled  to  a  predictive  sediment  diagenesis  model  that  computed 
sediment-water  fluxes  of  dissolved  oxygen,  methane,  ammonium,  nitrate,  and  phosphate,  based  on  com¬ 
puted  inputs  of  particulate  organic  matter.  Computed  values  of  deposition  and  fluxes  were  in  reasonable 
agreement  with  published  values.  Nutrient  loads  were  calculated  and  nutrient  budgets  were  constructed  as 
part  of  the  model  exercise.  Load  sources  included  river  inflows,  distributed  loads,  sewer  overflows,  and 
atmospheric  loading.  The  Sammamish  River  w'as  identified  as  the  largest  source  of  nutrients  to  Lake 
Washington,  followed  by  the  Cedar  River  and  other  distributed  sources.  The  majority  of  the  nutrient  load 
is  deposited  in  the  sediments.  A  lesser  amount  leaves  via  Lake  Union.  Nutrient  loads  in  this  study  were 
30  percent  (nitrogen)  to  60  percent  (phosphorus)  higher  than  the  loads  from  the  late  1970’s. 
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1 


Introduction 


Lake  Washington 

Lake  Washington  is  one  of  a  system  of  three  lakes  located  in  northwestern 
Washington  State  (Figure  1-1).  The  lake  provides  a  textbook  example  of 
eutrophication  remediation  through  reduction  of  nutrient  loads.  The  decline  and 
recovery  of  the  lake  are  the  subject  of  a  series  of  classic  works  by  W.  T. 
Edmondson,  his  students,  and  co-workers  (Edmondson  1972;  Edmondson  and 
Lehman  1981;  Edmondson  1991;  Edmondson  1994).  The  description  of  the  lake 
in  this  chapter  is  derived  from  these  sources  as  well  as  from  a  Web  site 
maintained  by  the  King  County  Department  of  Natural  Resources  and  Parks 
( http://dnr/metrokc.gov/wlr/waterres/lakes/Monitor.htm ). 

Lake  Washington  was  formed  at  the  end  of  an  ice  age  roughly  12,000  years 
ago.  Within  historic  times,  the  lake  had  one  major  inflow  at  its  northern  end,  the 
Sammamish  River,  which  drains  Lake  Sammamish.  The  outlet  was  the  Black 
River  at  the  southern  end.  The  Black  River  flowed  into  the  Cedar  River,  which 
led  to  Puget  Sound.  In  1916,  a  locked  ship  canal  was  constructed  that  connected 
Lake  Washington  with  Puget  Sound  through  Lake  Union  to  the  west.  The  level 
of  Lake  Washington  was  lowered  by  3.3  m,  and  the  Cedar  River  was  converted 
from  an  outflow  to  an  inflow.  The  major  outflow  was  now  through  the 
Chittenden  locks  at  the  western  end  of  Lake  Union. 

Lake  Washington  received  increasing  amounts  of  sewage  effluent  from  1941 
to  1963.  In  response,  concentrations  of  nutrients  and  phytoplankton  increased 
and  the  nature  of  the  phytoplankton  community  changed.  Between  1963  and 
1968,  sewage  discharge  to  the  lake  was  eliminated.  Concentrations  of  nutrients 
and  phytoplankton  and  the  proportion  of  blue-green  algae  all  declined  as  a  result. 
By  1975,  the  lake  had  recovered  from  its  eutrophic  state. 

Lake  Washington  is  the  largest  of  the  three  major  lakes  in  King  County,  and 
the  second  largest  natural  lake  in  the  state  of  Washington.  Lake  Washington’s 
two  major  influent  streams  are  the  Cedar  River  at  the  southern  end,  which 
contributes  about  57  percent  of  the  annual  runoff,  and  the  Sammamish  River  at 
the  northern  end,  which  contributes  27  percent  of  the  annual  runoff.  Total  annual 
flow  is  approximately  41  percent  of  the  lake  volume.  Residence  time  is  about 
2.4  years.  The  majority  of  the  immediate  watershed  (Figure  1-2)  is  highly 
developed  and  urban  in  nature  with  63  percent  fully  developed.  The  upper 
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portion  of  the  watershed  is  the  headwaters  of  the  Cedar  River  that  lie  in  the 
closed  Seattle  Water  Department  watershed. 
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The  basin  of  Lake  Washington  (Table  1-1)  is  a  deep,  narrow,  glacial  trough 
with  steeply  sloping  sides,  sculpted  by  the  Vashon  ice  sheet,  the  last  continental 
glacier  to  move  through  the  Seattle  area.  The  lake  is  6.3  m  (20.6  ft)  above  mean 
lower  low  tide  in  Puget  Sound,  to  which  it  is  connected  via  Lake  Union  and  the 
Lake  Washington  Ship  Canal.  Mercer  Island  lies  in  the  southern  half  of  the  lake, 
separated  from  the  east  shore  by  a  relatively  shallow  and  narrow  channel,  and 
from  the  west  shore  by  a  much  wider  and  deeper  channel. 


Table  1-1 

Physical  Characteristics  of  Lakes  Washington  and  Union 

Property 

Lake  Washington 

Lake  Union 

Drainage  area,  km2 

1,274 

1,554 

Surface  area,  km2 

87.6 

2.3 

Volume,  m3 

2  x  109 

2.5  x  107 

Mean  depth,  m 

32.9 

10 

Maximum  depth,  m 

65.2 

5 

Residence  time,  years 

2.3 

0.02 

Lake  Union 

Lake  Union  (Table  1-1)  was  formed  by  the  Vashon  glacier  about 
12,000  years  ago.  Originally  the  lake  was  separated  from  Lake  Washington  by 
the  land  extending  between  Portage  Bay  and  Union  Bay,  and  outflow  from  the 
lake  was  through  Salmon  Bay,  which  was  then  a  marine  inlet. 


Chapter  1  Introduction 


1-3 


In  1914,  the  hydrology  of  the  lake  was  significantly  altered  by  the 
construction  of  the  Fremont  and  Montlake  cuts  and  the  Hiram  M.  Chittenden 
Locks.  These  modifications  increased  inflow  to  Lake  Union  by  diverting  the 
outflow  from  Lake  Washington  through  Lake  Union  via  the  Montlake  Cut  and 
the  Ship  Canal.  Lake  Union  now  flushes  completely  about  once  a  week  at  high 
water  flows.  The  opening  of  the  Ship  Canal  allows  periodic  influx  of  salt  water 
from  Puget  Sound.  The  extent  of  the  saltwater  intrusion  varies  with  runoff 
volume  and  lock  operation.  During  the  rainy  season  and  spring  thaw,  runoff  from 
the  Cascade  foothills  is  high  and  the  lake  is  flushed  with  fresh  water.  As  the  flow 
drops  off  in  the  summer  and  boat  passage  through  the  lock  increases,  the 
intrusion  of  salt  water  through  the  locks  increases. 

The  intrusion  of  salt  water  forms  a  saline  wedge  that  flows  up  the  bottom  of 
the  Ship  Canal  and  the  lake.  The  saline  bottom  water  becomes  devoid  of  oxygen 
early  in  the  summer,  as  the  oxygen  is  used  by  bacteria  that  consume  organic 
matter.  The  anoxic  conditions  and  warm  water  temperatures  limit  the  amount  of 
habitat  available  to  fish.  The  irregular  influx  of  salt  water  poses  an  adaptation 
problem  for  the  benthic  organisms  that  form  the  basis  of  the  lake’s  food  web. 
Concern  for  these  effects  and  the  potential  for  permanent  damage  to  Lake 
Washington,  should  a  large  flow  of  salt  water  move  that  far,  prompted  the  U.S. 
Army  Corps  of  Engineers  in  1 964  to  build  a  saltwater  siphon  at  the  locks  to 
return  some  of  the  salt  water  to  Puget  Sound.  Beginning  in  the  early  1980’s  the 
Corps  of  Engineers  made  additional  modifications  in  lock  operation  to  help 
reduce  the  saline  intrusion  and  to  conserve  Lake  Washington  water  used  to  flush 
saline  water  from  the  locks. 

Lake  Union  is  unique  among  the  three  local  major  lakes  in  the  character  of  its 
watershed  (Figure  1-2).  It  is  the  most  heavily  urbanized  of  the  three,  draining 
residential,  commercial,  and  industrial  neighborhoods.  Its  shores  are  completely 
lined  by  marinas,  houseboat  moorage,  commercial  docks  and  dry  docks,  and 
industries.  In  the  past,  the  lake  received  sanitary  discharges  from  houseboats  and 
ships  and  industrial  discharges  from  businesses  along  the  shore,  as  well  as  fuel 
spills  and  discharges  from  ships  and  onshore  facilities.  While  pollution  inputs 
from  many  of  these  sources  have  decreased,  not  all  of  these  discharges  have 
stopped.  Combined-sewer  overflows  still  occur  and  are  generally  agreed  to  be  the 
worst  source  of  pollutants  and  pathogens  to  the  lake. 


Study  Objectives 

The  King  County  Water  Lands  Resources  Division  and  the  Wastewater 
Treatment  Division  are  in  the  process  of  developing  a  predictive  tool  to  assess 
water  quality  in  the  three  lakes,  the  Sammamish  River,  and  their  associated 
watersheds.  This  model  development  is  being  done  as  part  of  the  Sammamish- 
Washington  Analysis  and  Modeling  Program  (SWAMP).  The  assessment  has 
many  components,  but  generally  the  county  wants  to  develop  a  robust  tool  to 
assess  various  management  strategies  for  water  reuse  and  land  use.  Part  of  the 
assessment  tool  will  be  the  use  of  a  water  quality  computer  model.  The  SWAMP 
modeling  work  has  the  following  objectives: 
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a.  Understand  existing  conditions  in  the  water  bodies  of  the  watershed  and 
associated  risks  to  aquatic  life,  wildlife,  and  people. 

b.  Understand  future  conditions  in  the  water  bodies  of  the  watershed  and 
potential  risks  to  aquatic  life,  wildlife,  and  people  under  defined 
“buildout”  conditions. 

c.  Understand  the  effects  of  using  reclaimed  water  in  the  watershed  on 
existing  and  future  conditions  and  resulting  potential  risks  to  aquatic  life, 
wildlife,  and  people. 

In  September  2000,  the  King  County  Department  of  Natural  Resources  and 
Parks  signed  a  cooperative  agreement  with  the  U.S.  Army  Engineer  Research  and 
Development  Center  (ERDC)  to  provide  the  water  quality  models  necessary  for 
the  larger  SWAMP  effort.  ERDC  responsibilities  under  the  agreement  included 
the  following: 

a.  Configure  and  execute  a  model  for  Lake  Washington. 

b.  Submit  written  project  reports. 

c.  Provide  training  and  expertise  to  county  personnel. 

d.  Assist  in  assessing  model  operation,  interpretation  of  field  data,  and 
model  results. 

e.  Modify  code  as  required. 

f.  Provide  expertise  on  lake  and  river  dynamics  such  as  nutrient  dynamics, 
ecology,  and  sediment- water  interactions. 

For  the  modeling  effort,  ERDC  applied  the  CH3D-WES  hydrodynamic 
model  (Johnson  et  al.  1993)  and  the  CE-QUAL-ICM  eutrophication  model 
(Cerco  and  Cole  1993).  The  models  were  applied  to  Lakes  Washington  and 
Union  for  the  period  1995-1997  with  the  emphasis  on  Lake  Washington.  Lake 
Union  was  modeled  so  that  the  downstream  boundary  condition  could  be  located 
away  from  the  region  of  greatest  interest.  The  present  report  comprises  the 
primary  documentation  of  the  model  effort. 
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2  Flows,  Loads,  and 
Boundary  Conditions 


Flows 

Flows  into  Lake  Washington  may  be  divided  into  two  categories:  gauged  and 
distributed.  The  gauged  flows  enter  at  defined  locations  and  are  sufficiently  large 
to  merit  continuous  measurement  by  the  United  States  Geological  Survey 
(USGS).  Gauged  flows  (Figure  2-1)  include  the  Sammamish  River  near 
Woodinville  (12125200),  the  Cedar  River  at  Renton  (121 19000),  and  Thornton 
Creek  near  Seattle  (12128000).  Daily  gauged  flows  were  input  to  the  model  at 
corresponding  locations. 

The  distributed  flows  enter  at  undefined  locations  and  are  not  monitored.  The 
sponsor  provided  daily  distributed  flows  for  27  basins  (Figure  2-2).  These  were 
computed  based  on  precipitation  and  basin  characteristics.  The  distributed  flows 
were  pooled  into  10  locations  (Figure  2-1,  Table  2-1)  that  could  be  identified  as 
named  tributaries  or  features  and  input  to  the  model  on  a  daily  basis. 

Flows  into  Lake  Washington  showed  typical  seasonality.  Highest  flows 
occurred  from  November  through  March  or  April  (e.g..  Figure  2-3).  Lowest 
flows  occurred  in  August  or  September.  During  the  study  period,  flow  in  the 
Cedar  River  oscillated  around  the  long-term  average  in  1995  and  1996  and  was 
above  average  for  much  of  1997  (Figure  2-3).  Similar  behavior  was  observed  in 
the  Sammamish  River  (Figure  2-4).  Thornton  Creek  exhibited  above-average 
flows  for  the  1995-1996  wet  season  and  typical  flows  otherwise  (Figure  2-5). 

The  distributed  flows  showed  seasonality  consistent  with  the  gauged  flows 
(Figure  2-6).  No  distinct  wet  or  dry  years  were  apparent. 

For  the  study  period,  the  ranking  of  flow  sources  was  consistent  with 
statistics  reported  on  the  sponsor’s  Web  site  http://dm/metrokc.gov/wlr/waterres/ 
lakesMonitor.htm.  The  Cedar  River  provided  60  percent  of  the  runoff  followed 
by  the  Sammamish  River  with  28  percent.  Thornton  Creek  provided 
approximately  1  percent  and  was  revealed  to  be  a  minor  source  despite  the 
presence  of  a  USGS  gauge.  Distributed  flows  contributed  the  remaining 
1 1  percent. 
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Loads  of  nutrients,  organic  matter,  suspended  solids,  and  pathogens  enter 
Lake  Washington  from  three  primary  sources:  runoff,  combined-sewer  overflows 
(CSOs),  and  atmospheric  deposition.  Point-source  loads  from  municipal  waste 
were  eliminated  in  the  period  from  1963  to  1968. 

Runoff 

Daily  runoff  loads  were  computed  as  the  product  of  flow  and  concentration. 
Concentrations  in  gauged  flows,  observed  at  approximately  monthly  intervals, 
were  provided  by  the  sponsor.  Concentrations  were  reported  for  the  following: 

•  Total  Nitrogen 

•  Total  Phosphorus 

•  Total  Suspended  Solids 

•  Nitrate+Nitrite  Nitrogen 

•  Ortho-Phosphorus 

•  Ammonium 

•  Fecal  Coliform 

Linear  interpolation  was  used  to  create  daily  concentrations  from  the 
monthly  values.  The  record  for  the  Cedar  River  commenced  in  April  1996.  Prior 
to  this  date,  monthly-average  values  from  the  period  of  record  were  employed. 
Daily  concentrations  were  multiplied  by  daily  flows  to  create  a  time  series  of 
daily  loads. 

Nitrogen  and  phosphorus  concentrations  in  distributed  flows  from  each  basin 
were  provided  by  the  sponsor.  These  were  obtained  from  regression  relationships 
based  on  precipitation  and  impervious  area.  Daily  concentrations  were  multiplied 
by  daily  flows  to  create  a  time  series  of  daily  loads.  Loads  were  input  at  one 
location  for  each  of  the  basins.  Locations  corresponded  with  locations  of 
distributed  flows  when  these  were  present. 

Fecal  coliform  concentrations  in  distributed  flows  were  based  on 
observations  provided  by  the  sponsor.  Constant  values,  ~300  mpn/100  mL 
(where  mpn  =  most  probably  number),  were  employed  for  each  basin.  Constant 
suspended  solids  concentrations,  5  to  9  g  m'3,  were  similarly  derived  and 
employed.  Consequently,  variation  in  daily  loads  from  each  basin  reflects 
variation  in  flow  only. 

Organic  carbon.  A  limited  quantity  of  organic  carbon  observations, 
collected  in  2002,  was  available  for  the  Cedar  and  Sammamish  Rivers.  Mean 
total  organic  carbon  concentrations  were  6.2  g  m'3  for  the  Sammamish  and 
2.4  g  m'3  for  the  Cedar. 

This  situation  of  limited  organic  carbon  observations  is  common. 
Consequently,  organic  carbon  loads  for  the  model  are  computed  as  a  multiple  of 
organic  nitrogen  loads.  Based  on  the  observations,  a  carbon-to-nitrogen  ratio  of 
20  was  used. 
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Mapping  to  model  variables.  The  model  requires  loads  split  into  multiple 
forms  of  organic  carbon,  nitrogen,  and  phosphorus  (Table  2-2).  Organic  nitrogen 
and  phosphorus  were  estimated  as  the  difference  between  reported  total  and 
mineral  forms.  The  dissolved  fraction  of  organic  carbon  and  nitrogen  was  derived 
from  the  organic  carbon  observations.  The  dissolved  fraction  of  organic 
phosphorus  was  derived  from  data  reported  by  Edmondson  and  Lehman  (1981). 

Combined-sewer  overflows 

Combined-sewer  overflows,  sanitary' -sewer  overflows,  and  spills  were 
reported  at  28  locations  (Figure  2-7)  during  the  study  period.  Daily  flow  volume 
for  each  incident  was  provided  by  the  sponsor.  The  sponsor  also  provided 
characteristic  concentrations  for  the  discharges  (Table  2-3).  Concentrations  in 
CSOs  were  based  on  median  concentrations  observed  in  the  Landers  CSO. 
Concentrations  in  the  sanitary-sewer  overflows  and  spills  were  based  on  median 
concentrations  in  the  West  Point  Sewage  Treatment  Plant  influent.  Daily  loads 
were  computed  as  the  product  of  flow  and  concentration  and  input  to  the  model 
at  appropriate  locations.  All  organic  substances  were  split  25  percent  dissolved, 
50  percent  labile  particulate,  and  25  percent  refractory  particulate.  Flows  were 
considered  negligible  relative  to  runoff  volume  and  were  not  included  in  the 
model. 

Atmospheric  deposition 

Atmospheric  loads  were  obtained  from  Edmondson  and  Lehman  (1981). 
Their  mean  annual  loads  were  converted  to  daily  areal  loads  for  input  to  the 
model.  Total  nitrogen  loading  of  1.08  mg  m'2  day1  was  split  into  50  percent 
nitrate,  46  percent  dissolved  organic  nitrogen,  and  4  percent  ammonium.  Total 
phosphorus  loading  of  0.27  mg  m'2  day'1  was  split  into  equal  portions  of 
phosphate  and  dissolved  organic  phosphorus. 

Load  summary 

Nitrogen  and  phosphorus.  Time  series  of  nitrogen  (Figures  2-8  through 
2-1 1)  and  phosphorus  (Figures  2-13  through  2-16)  loads  corresponded  closely 
with  time  series  of  flow  (Figures  2-3  through  2-6).  The  correspondence  of  flow 
and  load  indicated  little  dependence  of  concentration  on  flow.  The  general  lack 
of  correspondence  between  flow  and  concentration  was  previously  noted  by 
Edmondson  and  Lehman  (1981).  Nitrate  was  the  predominant  component  of  the 
nitrogen  load  while  phosphate  composed  less  than  half  the  phosphorus  load. 

The  ranking  of  loads  from  various  runoff  sources  did  not  follow  the  ranking 
of  flows.  The  Sammamish  River  was  the  predominant  nutrient  source  (Table  2-7) 
even  though  flow  was  higher  in  the  Cedar.  The  proportion  of  nutrient  loads  from 
distributed  sources,  18  percent  to  25  percent,  exceeded  the  proportion  of  runoff 
from  distributed  sources,  1 1  percent.  The  Sammamish  contributed  greater  loads 
because  this  river  was  nutrient-enriched  relative  to  the  Cedar  (Tables  2-4,  2-5). 
Nutrient  concentrations  in  distributed  loads,  as  exemplified  by  Thornton  Creek 
(Table  2-6),  were  enriched  relative  to  both  large  rivers. 
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An  enormous  CSO  event  occurred  at  the  end  of  1996  and  extended  into  1997 
(Figures  2-12,  2-17).  Loads  from  this  event  were  comparable  to  the  highest  daily 
loading  from  major  runoff  sources.  The  event  was  unique  to  the  study  period. 
Other  overflows  and  spills  were  of  much  less  magnitude  and,  on  an  annual  basis, 
loads  from  these  sources  were  minor  relative  to  major  runoff  sources  (Tables  2-7, 
2-8). 


Atmospheric  deposition  was  a  minor  component  of  the  annual  total  nitrogen 
budget.  Atmospheric  phosphorus  deposition  was  less  than  major  runoff  sources 
but  comprised  a  potentially  significant  load,  especially  during  periods  of  low 
runoff. 

Fecal  coliform.  Fecal  coliform  loads  in  runoff  (Figures  2-18  through  2-22) 
showed  less  temporal  correspondence  to  flows  than  did  nutrient  loads.  Loads  in 
the  Cedar  were  “spiky”  with  multiple  peaks  throughout  the  study  period 
(Figure  2-18)  while  loads  in  the  Sammamish  were  dominated  by  two  events  that 
occurred  in  late  1995  and  in  late  1996  (Figure  2-19).  On  an  annual  basis,  the 
major  loads  were  from  the  Sammamish  and  Cedar,  followed  by  distributed  loads 
(Table  2-9).  CSO  loads  included  several  events  that  rivaled  loads  from  major 
runoff  sources  (Figure  2-22),  but  spills  and  overflows  were  minor  sources  on  a 
long-term  basis  (Table  2-9). 

Suspended  solids.  As  with  nutrients,  the  time  series  of  suspended  solids  in 
major  runoff  sources  corresponded  with  flows  (Figures  2-23  through  2-27). 
Suspended  solids  concentrations  in  the  Cedar  (Table  2-4)  exceeded  the 
Sammamish  (Table  2-5)  so  that  the  Cedar  was  the  dominant  source  to  the  system 
(Table  2-10).  Distributed  sources,  exemplified  by  Thornton  Creek  (Table  2-6), 
were  enriched  with  solids  relative  to  major  runoff  sources;  but  the  role  of 
distributed  solids  sources  (Table  2-10)  was  less  than  distributed  nutrient  sources. 
The  late  1996-early  1997  CSO  event  (Figure  2-27)  produced  solids  loads  that 
rivaled  the  major  rivers;  but  overflows  and  spills  were  a  minor  component  of  the 
long-term  suspended  solids  budget  (Table  2-10). 


Reference 

Edmondson,  W.  T.,  and  Lehman,  J.  (1981).  “The  effect  of  changes  in  the  nutrient 
income  on  the  condition  of  Lake  Washington,”  Limnology  and 
Oceanography  26(1),  1-29. 
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Table  2-1 

Flows  and  Loads 

Source 

Type 

Inflow  Location 

Load 

Sammamish  River 

Gauged 

Yes 

Observed  flow  *  interpolated 
concentration 

Cedar  River 

Yes 

Thornton  Creek 

Yes 

Lyons 

Distributed 

Yes,  combined  with 
McAleer 

Computed  flow  *  computed 
concentration 

McAleer 

Lake  Park  Shore 

No 

Magnuson  Park 

No 

Ravenna 

No 

No 

Madison  Park 

No 

Ranier-Genesee 

No 

Seward  Park 

No 

Lakeridge 

No 

Finn-Ken 

Yes,  combined  with 

Juanita 

Juanita 

Forbes 

Yes,  combined  with 

Kirkland 

Kirkland  Shores 

Yarrow  Bay 

Yes,  combined  with  Cozy 
Cove 

Cozy  Cove 

Fairweather 

Yes,  combined  with 
Meydenbaur 

South  of 

Evergreen  Pt. 

No 

Kelsey-Mercer 

Yes 

Chism  Beach 

No 

West  Mercer 

Yes 

Coal  Creek 

Yes 

May  Creek 

Yes,  combined  with 
Pleasure  Pt. 

Pleasure  Pt. 

Coulon  Park 

Yes,  combined  with 

Cedar  River 

East  Mercer  Island 

No 

Table  2-2 

Sources  and  Splits  of  Organic  Matter 

Obtained  from 

Split  into 

Carbon 

Organic  nitrogen  *  20 

85  percent  dissolved 

7.5  percent  labile  particulate 

7.5  percent  refractory  particulate 

■ 

Total  nitrogen  minus  nitrate  minus  ammonium 

85  percent  dissolved 

7.5  percent  labile  particulate 

7.5  percent  refractory  particulate 

Phosphorus 

Total  phosphorus  minus  phosphate 

10  percent  dissolved 

45  percent  labile  particulate 

45  percent  refractory  particulate 
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Table  2-3 

Concentrations  in  Sewage  Overflows 

Type  of 
Overflow 

no23> 
g  m"3 

nh4, 

-3 

g  rn 

Total  N, 

„  ^-3 
g  m 

Phosphate, 

-3 

gm 

Total  P, 
g  m"3 

Total 

Organic 

Carbon, 

-3 

g  m 

Total 

Suspended 
Solids  (TSS), 
g  rrf3 

Fecal 
Coliform, 
mpn/100  mL 

Combined 

Sewers 

0.51 

1.1 

2.58 

0.1 

0.42 

10.2 

22.9 

24,000 

Sanitary 

Sewers 

0.51 

8.95 

20.7 

2.22 

2.85 

40 

130 

24,000 

Table  2-4 

Cedar  River  Concentrations 

Total  N, 

no3, 

nh4, 

Total  P, 

Fecal  Coliform, 

HSi 

Concentration 

-3 

g  m 

g  m"3 

-3 

g  rn 

~3 

g  m 

IH 

mpn/100  mL 

1™ 

Mean 

0.42 

0.27 

0.030 

0.036 

0.011 

255 

11.9 

Median 

0.38 

0.25 

0.028 

0.024 

0.008 

110 

5.1 

Max 

1.12 

0.62 

0.071 

0.199 

0.061 

2900 

57.2 

Min 

0.20 

0.13 

0.010 

0.008 

0.003 

2 

1.5 

Stdev 

0.20 

0.10 

0.014 

0.033 

0.011 

514 

14.5 

N 

41 

41 

30 

40 

41 

39 

41 

Table  2-5 

Sammamish  River  Concentrations 

Concentration 

Total  N, 

-3 

g  rn 

no3, 

-3 

g  m 

nh4, 

-3 

g  m 

Total  P, 

-3 

g  rn 

P04, 

-3 

gm 

Fecal  Coliform, 
mpn/100  mL 

TSS, 

-3 

gm 

Mean 

0.86 

0.52 

0.042 

0.057 

0.023 

286 

6.3 

Median 

0.81 

0.45 

0.038 

0.054 

0.022 

205 

4.9 

Max 

1.49 

0.97 

0.092 

0.157 

0.044 

3200 

44.9 

Min 

0.57 

0.21 

0.018 

0.032 

0.010 

27 

wmi 

Stdev 

0.22 

0.18 

0.019 

0.022 

0.009 

468 

N 

48 

48 

44 

48 

48 

50 

48  || 

2-6 


Chapter  2  Flows,  Loads,  and  Boundary  Conditions 


Table  2-6 

Thornton  Creek  Concentrations 

Concentration 

mm 

no3, 

-3 

g  rn 

nh4, 

g  m'3 

Total  P, 

„  -3 

g  m 

P04, 

-3 

g  rn 

Fecal  Coliform, 
mpn/100  mL 

TSS, 
g  m'3 

Mean 

1.46 

1.14 

0.048 

0.106 

0.040 

1782 

23.4 

Median 

1.50 

1.29 

0.048 

0.080 

0.035 

780 

7.0 

Max 

2.28 

2.03 

0.113 

0.340 

0.101 

7900 

170.0 

Min 

0.82 

0.36 

0.020 

0.007 

0.014 

150 

1.3 

Stdev 

0.30 

0.41 

0.023 

0.074 

0.021 

2011 

38.1 

N 

48 

48 

45 

48 

48 

48 

48 

Table  2-7 

Nitrogen  Loads,  kg  day'1 

Source 

1995 

1996 

1997 

Fraction 

Cedar 

709 

841 

1031 

28.4 

Sammamish 

1092 

1354 

1308 

41.3 

Thornton 

55 

66 

44 

1.8 

Distributed 

614 

894 

760 

24.9 

cso 

2 

33 

11 

0.5 

Atmospheric 

95 

95 

95 

3.1 

Table  2-8 

Phosphorus  Loads,  kg  day'1 

Source 

1995 

1996 

1997 

Fraction 

Cedar 

42.0 

62.5 

85.9 

30.5 

Sammamish 

82.1 

76.7 

68.8 

36.5 

Thornton 

4.0 

3.0 

3.1 

1.6 

Distributed 

33.0 

43.6 

38.3 

18.4 

CSO 

0.4 

4.8 

1.7 

1.1 

Atmospheric 

24.7 

24.7 

24.7 

11.9 

Table  2-9 

Fecal  Coliform  Loads,  10s  organisms  day'1 

Source 

1995 

1996 

1997 

Fraction 

Cedar 

1 ,820,638 

1,960,678 

5,110,896 

31 .2 

bbbhh 

4,662,155 

5,538,131 

1 ,589,455 

41 .4 

579,672 

356,444 

475,251 

5.0 

Distributed 

1,229,815 

1 ,689,904 

1,501,836 

15.5 

CSO 

218,180 

1,202,159 

548,756 

6.9 
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Table  2-10 

Suspended  Solids  Loads,  kg  day'1 


Source 

1995 

1996 

1997 

Fraction 

Cedar 

12684 

19520 

23553 

53.7 

Sammamish 

15484 

11717 

9315 

35.2 

Thornton 

1174 

484 

558 

2.1 

Distributed 

2484 

3443 

3038 

8.6 

CSO 

21 

235 

83 

0.3 

Figure  2-1.  Locations  of  gauged  and  distributed  flows 
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Figure  2-2.  Distributed  flow  basins 
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jre  2-3.  Flow  in  Cedar  River 
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Figure  2-4.  Flow  in  Sammamish  River 
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cubic  meters  /  second 


Figure  2-6.  Distributed  flows 
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Figure  2-7.  Location  of  sewer  overflows  and  spills 
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Figure  2-19.  Daily  fecal  coliform  load  from  Sammamish  River 
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Figure  2-23.  Daily  suspended  solids  load  from  Cedar  River 


2-20 


Chapter  2  Flows,  Loads,  and  Boundary  Conditions 


Figure  2-24.  Daily  suspended  solids  load  from  Sammamish  River 
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Figure  2-25.  Daily  suspended  solids  load  from  Thornton  Creek 
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Figure  2-26.  Daily  suspended  solids  load  from  distributed  sources 


Figure  2-27.  Daily  suspended  solids  load  from  sewer  overflows  and  spills 
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3  Water  Quality  Model 
Formulation 


Introduction 

CE-QUAL-ICM  was  designed  to  be  a  flexible,  widely  applicable 
eutrophication  model.  Initial  application  was  to  Chesapeake  Bay  (Cerco  and  Cole 
1994).  Subsequent  additional  applications  included  the  Delaware  inland  bays 
(Cerco  et  al.  1994),  Newark  Bay  (Cerco  and  Bunch  1997),  and  the  San  Juan  Bay 
Estuary  (Bunch  et  al.  2000).  Each  model  application  employed  a  different 
combination  of  model  features  and  required  addition  of  system-specific 
capabilities.  This  chapter  describes  general  features  and  site-specific 
developments  of  the  model  as  presently  applied  to  the  water  column  of  Lake 
Washington. 


Conservation  of  Mass  Equation 

The  foundation  of  CE-QUAL-ICM  is  the  solution  to  the  three-dimensional 
mass-conservation  equation  for  a  control  volume.  Control  volumes  correspond  to 
cells  on  the  model  grid.  CE-QUAL-ICM  solves,  for  each  volume  and  for  each 
state  variable,  the  equation: 

S  v  •  C  ”  »  SC 

Vt  ^-=1  Qk  •  C*+I A*  •  Dk  •  +  (3'1} 

8  t  *-j  k=i  8  xk 


in  which 


Vj  =  volume  of fl  control  volume  (m3) 

Cj  =  concentration  in  fl  control  volume  (g  m'3) 
t,x-  temporal  and  spatial  coordinates 
n  =  number  of  flow  faces  attached  to  f‘  control  volume 
Qk  =  volumetric  flow  across  flow  face  k  of f'  control  volume  (m3  s'1) 
Ck  =  concentration  in  flow  across  face  k  (g  m'3) 

Ak  =  area  of  flow  face  k  (m2) 
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Die  =  diffusion  coefficient  at  flow  face  k  (m2  s'1) 

Si  =  external  loads  and  kinetic  sources  and  sinks  in  fh  control  volume 
(gs'1) 

Solution  of  Equation  3-1  on  a  digital  computer  requires  discretization  of  the 
continuous  derivatives  and  specification  of  parameter  values.  The  equation  is 
solved  using  the  QUICKEST  algorithm  (Leonard  1979)  in  the  horizontal  plane 
and  a  Crank -Nicolson  scheme  in  the  vertical  direction.  Discrete  time-steps, 
determined  by  computational  stability  requirements,  are  approximately 
15  minutes. 


State  Variables 

As  applied  to  Lake  Washington,  the  CE-QUAL-ICM  model  incorporates  20 
state  variables  in  the  water  column  including  physical  variables,  phytoplankton, 
multiple  forms  of  carbon,  nitrogen,  and  phosphorus,  a  pathogen,  and  two 
toxicants  (Table  3-1).  The  pathogen  and  toxicants  were  added  specifically  for  the 
Lake  Washington  application  and  are  described  in  a  subsequent  chapter. 


Table  3-1 

Lake  Washington  Water  Quality  Model  State  Variables 

Temperature 

Fixed  Solids 

Phytoplankton 

Fecal  Coliform 

Dissolved  Organic  Carbon 

Labile  Particulate  Organic  Carbon 

Refractory  Particulate  Organic  Carbon 

Ammonium 

Nitrate+Nitrite 

Dissolved  Organic  Nitrogen 

Labile  Particulate  Organic  Nitrogen 

Refractory  Particulate  Organic  Nitrogen 

Total  Phosphate 

Dissolved  Organic  Phosphorus 

Labile  Particulate  Organic  Phosphorus 

Refractory  Particulate  Organic  Phosphorus 

Chemical  Oxygen  Demand 

Dissolved  Oxygen 

Toxicant  One 

Toxicant  Two 

Algae 

The  coding  of  CE-QUAL-ICM  allows  the  specification  of  up  to  three  algal 
groups.  One  group  is  presently  activated  for  Lake  Washington. 

Zooplankton 

The  coding  of  CE-QUAL-ICM  allows  the  specification  of  two  zooplankton 
groups:  microzooplankton  and  mesozooplankton.  These  are  not  presently 
activated. 

Organic  carbon 

Three  organic  carbon  state  variables  are  considered:  dissolved,  labile 
particulate,  and  refractor}'  particulate.  Labile  and  refractor}'  distinctions  are  based 
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upon  the  time  scale  of  decomposition.  Labile  organic  carbon  decomposes  on  a 
time  scale  of  days  to  weeks  while  refractory  organic  carbon  requires  more  time. 
Labile  organic  carbon  decomposes  rapidly  in  the  water  column  or  the  sediments. 
Refractory  organic  carbon  decomposes  slowly,  primarily  in  the  sediments,  and 
may  contribute  to  sediment  oxygen  demand  years  after  deposition. 

Nitrogen 

Nitrogen  is  first  divided  into  available  and  unavailable  fractions.  Available 
refers  to  employment  in  algal  nutrition.  Two  available  forms  are  considered: 
reduced  and  oxidized  nitrogen.  For  the  Lake  Washington  application,  reduced 
nitrogen  consists  exclusively  of  ammonium.  Nitrate  and  nitrite  compose  the 
oxidized  nitrogen  pool.  Both  reduced  and  oxidized  nitrogen  are  utilized  to  fulfill 
algal  nutrient  requirements.  The  primary  reason  for  distinguishing  the  two  is  that 
ammonium  is  oxidized  by  nitrifying  bacteria  into  nitrite  and,  subsequently, 
nitrate.  This  oxidation  can  be  a  significant  sink  of  oxygen  in  the  water  column 
and  sediments. 

Unavailable  nitrogen  state  variables  are  dissolved  organic  nitrogen,  labile 
particulate  organic  nitrogen,  and  refractory  particulate  organic  nitrogen. 

Phosphorus 

As  with  nitrogen,  phosphorus  is  first  divided  into  available  and  unavailable 
fractions.  Only  a  single  available  form,  dissolved  phosphate,  is  considered.  The 
model  framework  allows  for  exchange  of  phosphate  between  dissolved  and 
particulate  (sorbed  to  solids)  forms,  but  this  option  is  not  implemented  in  the 
present  application.  Three  forms  of  unavailable  phosphorus  are  considered: 
dissolved  organic  phosphorus,  labile  particulate  organic  phosphorus,  and 
refractory  particulate  organic  phosphoms. 

Silica 

Silica  is  included  in  the  complete  version  of  CE-QUAL-ICM  but  is  not 
activated  in  the  Lake  Washington  application. 

Chemical  oxygen  demand 

Chemical  oxygen  demand  is  the  concentration  of  reduced  substances  that  are 
oxidized  by  abiotic  processes.  In  fresh  water,  a  primary  component  of  chemical 
oxygen  demand  is  methane  released  from  sediments.  Oxidation  of  methane  may 
remove  substantial  quantities  of  dissolved  oxygen  from  the  water  column. 

Dissolved  oxygen 

Dissolved  oxygen  is  required  for  the  existence  of  higher  life  forms.  Oxygen 
availability  determines  the  distribution  of  organisms  and  the  flows  of  energy  and 
nutrients  in  an  ecosystem.  Dissolved  oxygen  is  a  central  component  of  the  water- 
quality  model. 
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Salinity 

Salinity  is  included  in  the  coding  of  CE-QUAL-ICM  but  is  not  activated  in 
Lake  Washington  at  present. 

Temperature 

Temperature  is  a  primary  determinant  of  the  rate  of  biochemical  reactions. 
Reaction  rates  increase  as  a  function  of  temperature  although  extreme 
temperatures  may  result  in  the  mortality  of  organisms  and  a  decrease  in  kinetics 
rates. 

Fixed  solids 

Fixed  solids  are  the  mineral  fraction  of  total  suspended  solids.  Solids  are 
considered  primarily  for  their  role  in  light  attenuation. 

The  remainder  of  this  chapter  is  devoted  to  detailing  the  kinetics  sources  and 
sinks  and  to  reporting  parameter  values.  For  notational  simplicity,  the  transport 
terms  are  dropped  in  the  reporting  of  kinetics  formulations. 


Algae 

Algal  sources  and  sinks  in  the  conservation  equation  include  production, 
metabolism,  predation,  and  settling.  These  are  expressed  as 


f 


8  t 


B  = 


G-BM-Wa  •  — 
dz 


B-PR 


in  which 

B  =  algal  biomass,  expressed  as  carbon  (g  C  m'3) 
G  =  growth  (day-1) 

BM=  basal  metabolism  (day-1) 

Wa  =  algal  settling  velocity  (m  day'1) 

PR  =  predation  (g  C  m'3  day'1) 
z  =  vertical  coordinate 


(3-2) 


Production 

Production  by  phytoplankton  is  determined  by  the  intensity  of  light,  by  the 
availability  of  nutrients,  and  by  the  ambient  temperature. 

Light 

The  influence  of  light  on  phytoplankton  production  is  represented  by  a 
chlorophyll-specific  production  equation  (Jassby  and  Platt  1976): 
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pB  =  pBm  —r - 

V/2+/£2 


in  which: 

Pn  =  photosvnthetic  rate  (g  C  g'1  Chi  day1) 

P,sm  =  maximum  photosynthetic  rate  (g  C  g'1  Chi  day'1) 

1  =  irradiance  (E  m'2  day'1) 

Parameter  Ik  is  defined  as  the  irradiance  at  which  the  initial  slope  of  the 
production  versus  irradiance  relationship  (Figure  3-1)  intersects  the  value  of  F^m: 


L *  =  ^ 
a 


(3-4) 


in  which  a  is  the  initial  slope  of  production  versus  irradiance  relationship  (g  C  g'1 
Chi  (E  m'2)'1). 


Figure  3-1 .  Production  versus  irradiance  curve 


Chlorophyll-specific  production  rate  is  readily  converted  to  carbon-specific 
growth  rate,  for  use  in  Equation  3-2,  through  division  by  the  carbon-to- 
chlorophyll  ratio: 
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G=  p^_ 

CChl 


(3-5) 


in  which  CChl  is  the  carbon-to-chlorophyll  ratio  (g  C  g'1  chlorophyll  a). 

Nutrients 

Carbon,  nitrogen,  and  phosphorus  are  the  primary  nutrients  required  for  algal 
growth.  Inorganic  carbon  is  usually  available  in  excess  and  is  not  considered  in 
the  model.  The  effects  of  the  remaining  nutrients  on  growth  are  described  by  the 
formulation  commonly  referred  to  as  Michaelis-Menton  kinetics  (Figure  3-2): 


M- 


D 

KHd  +  D 


(3-6) 


in  wrhich 

f(N)  -  nutrient  limitation  on  algal  production  (0  <f(N)  <  1) 
D  =  concentration  of  dissolved  nutrient  (g  m'3) 

KHd  =  half-saturation  constant  for  nutrient  uptake  (g  m'3) 


_ N/KH _ 

Figure  3-2.  Michaelis-Menton  formulation  for  nutrient-limited  growth 

Temperature 

Algal  production  increases  as  a  function  of  temperature  until  an  optimum 
temperature  or  temperature  range  is  reached.  Above  the  optimum,  production 
declines  until  a  temperature  lethal  to  the  organisms  is  attained.  Numerous 
functional  representations  of  temperature  effects  are  available.  Inspection  of 
growth  versus  temperature  data  indicates  a  function  similar  to  a  Gaussian 
probability  curve  (Figure  3-3)  provides  a  good  fit  to  observations: 
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f(T)  =  eKTg!  *  a -Toptf  when  T  <  Topt 

=  e-KTg2.(Topt-TfwhenT  >  Top{ 


(3-7) 


in  which 

T  =  temperature  (°C) 

Topt  =  optimal  temperature  for  algal  growth  (°C) 

KTgl  =  effect  of  temperature  below  Topt  on  growth  (°C'2) 
KTg2  =  effect  of  temperature  above  Topt  on  growth  (°C'2) 


Figure  3-3.  Relation  of  algal  production  to  temperature 

Constructing  the  photosynthesis  versus  irradiance  curve 

A  production  versus  irradiance  relationship  is  constructed  for  each  model  cell 
at  each  time-step.  First,  the  maximum  photosynthetic  rate  under  ambient 
temperature  and  nutrient  concentrations  is  determined: 

PB  m(N,T)  =  pB  m  •  f(T)  •  (3-8) 

in  which  F^mfN,  T)  is  the  maximum  photosynthetic  rate  under  ambient 
temperature  and  nutrient  concentrations  (g  C  g'1  Chi  day'1).  The  single  most 
limiting  nutrient  is  employed  in  determining  the  nutrient  limitation. 

Next,  parameter  Ik  is  derived  from  Equation  3-4.  Finally,  the  production 
versus  irradiance  relationship  is  constructed  using  T)  and  Ik.  The  resulting 

production  versus  irradiance  curve  exhibits  three  regions  (Figure  3-4).  For  /» 

Ik ,  the  value  of  the  term  1  /  (f  +  Ik2)  *  approaches  unity  and  temperature  and 
nutrients  are  the  primary  factors  that  influence  production.  For  1 «  Ik , 
production  is  determined  solely  by  a  and  irradiance  I.  In  the  region  where  the 


Chapter  3  Water  Quality  Model  Formulation 


3-7 


initial  slope  of  the  production  versus  irradiance  curve  intercepts  the  line 
indicating  production  at  optimal  illumination,  1  ~  Ik,  production  is  determined  by 
the  combined  effects  of  temperature,  nutrients,  and  light. 


Figure  3-4.  Effects  of  light  and  nutrients  on  production  versus  irradiance  curve, 
calculated  for  a  =  8  (g  C  g"1  Chi  (E  m'2)'1) 

Irradiance 

Solar  radiation  at  the  water  surface  is  input  to  the  model  as  part  of  the  heat 
transfer  calculations.  An  empirical  conversion  is  employed  to  convert  between 
short-wave  solar  radiation  and  photosynthetically  active  radiation  (Einstein’s  m2 
=  0.143  •  watts  m'2). 

Respiration 

Two  forms  of  respiration  are  considered  in  the  model:  photorespiration  and 
basal  metabolism.  Photorespiration  represents  the  energy  expended  by  carbon 
fixation  and  is  a  fixed  fraction  of  production.  In  the  event  of  no  production  (e.g., 
at  night),  photorespiration  is  zero.  Basal  metabolism  is  a  continuous  energy 
expenditure  to  maintain  basic  life  processes.  In  the  model,  metabolism  is 
considered  to  be  an  exponentially  increasing  function  of  temperature 
(Figure  3-5).  Total  respiration  is  represented: 


R  =  Presp  •  G  +  BM  •  eKTb  *  <T'Tr^) 


(3-9) 


in  which 


3-8 


Chapter  3  Water  Quality  Model  Formulation 


Presp  -  photorespiration  (0  <Presp  <  1) 

BM=  metabolic  rate  at  reference  temperature  Tref  (day'1) 
KTb  =  effect  of  temperature  on  metabolism  (°C1) 

Tref  =  reference  temperature  for  metabolism  (°C) 


T  -  Tref 


Figure  3-5.  Exponential  temperature  relationship  employed  for  metabolism  and 
other  processes 


Predation 

The  predation  term  includes  the  activity  of  zooplankton,  filter-feeding 
benthos,  and  other  pelagic  filter  feeders  including  planktivorous  fish.  A  variety  of 
formulations  are  available  to  represent  predation.  For  Lake  Washington,  a 
piecewise  linear  predation  term  worked  well: 

PR  =  BPR»B  (3-10) 

in  which  BPR  is  the  base  predation  rate  (day"1). 

Accounting  for  algal  phosphorus 

The  amount  of  phosphorus  incorporated  in  algal  biomass  is  quantified 
through  a  stoichiometric  ratio.  Thus,  total  phosphoms  in  the  model  is  expressed: 

TotP  =  P04  +  APC  »B  +  DOP  +  LPOP  +  RPOP  (3-11) 

in  which 

TotP  =  total  phosphorus  (g  P  m"3) 

P04  =  dissolved  phosphate  (g  P  m"3) 


Chapter  3  Water  Quality  Model  Formulation 


3-9 


APC  =  algal  phosphorus-to-carbon  ratio  (g  P  g'1  C) 

DOP  =  dissolved  organic  phosphorus  (g  P  m'3) 

LPOP  =  labile  particulate  organic  phosphorus  (g  P  m'3) 

RPOP  -  refractory  particulate  organic  phosphorus  (g  P  m'3) 

Algae  take  up  dissolved  phosphate  during  production  and  release  dissolved 
phosphate  and  organic  phosphorus  through  respiration.  The  fate  of  phosphorus 
released  by  respiration  is  determined  by  empirical  distribution  coefficients.  The 
fate  of  algal  phosphorus  recycled  by  predation  is  determined  by  a  second  set  of 
distribution  parameters. 

Accounting  for  algal  nitrogen 

Model  nitrogen  state  variables  include  ammonium,  nitrate+nitrite,  dissolved 
organic  nitrogen,  labile  particulate  organic  nitrogen,  and  refractory  particulate 
organic  nitrogen.  The  amount  of  nitrogen  incorporated  in  algal  biomass  is 
quantified  through  a  stoichiometric  ratio.  Thus,  total  nitrogen  in  the  model  is 
expressed: 

TotN  =  NH4  +  N023  +  ANC  •  B  +  DON + LPON  +  RPON  (3-12) 

in  which: 

TotN  =  total  nitrogen  (g  N  m'3) 

NH4  =  ammonium  (g  N  m'3) 

NO23  =  nitrate+nitrite  (g  N  m'3) 

ANC  =  algal  nitrogen-to-carbon  ratio  (g  N  g'1  C) 

DON  =  dissolved  organic  nitrogen  (g  N  m'3) 

LPON  =  labile  particulate  organic  nitrogen  (g  N  m'3) 

RPON  =  refractory  particulate  organic  nitrogen  (g  N  m'3) 

As  with  phosphorus,  the  fate  of  algal  nitrogen  released  by  metabolism  and 
predation  is  represented  by  distribution  coefficients. 

Algal  nitrogen  preference 

Algae  take  up  ammonium  and  nitrate+nitrite  during  production  and  release 
ammonium  and  organic  nitrogen  through  respiration.  Nitrate+nitrite  is  internally 
reduced  to  ammonium  before  synthesis  into  biomass  occurs  (Parsons  et  al.  1984). 
Trace  concentrations  of  ammonium  inhibit  nitrate  reduction  so  that,  in  the 
presence  of  multiple  nitrogenous  nutrients,  ammonium  is  utilized  first.  The 
“preference”  of  algae  for  ammonium  is  expressed  by  an  empirical  function 
(Thomann  and  Fitzpatrick  1982): 
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PN  =  NH4* 


_ NO  23 _ 

(KHn  +  NH4  )  •  (KHn  +  N023 ; 


+nh4. 


_ KHn _ 

(NH4  +  N023 )  •  (KHn  +  N023 ) 


in  which 

PN  =  algal  preference  for  ammonium  uptake  (0  <Pn  <  1) 

KHn  =  half-saturation  concentration  for  algal  nitrogen  uptake  (g  N  m'3) 

The  function  has  two  limiting  values  (Figure  3-6).  When  nitrate+nitrite  is 
absent,  the  preference  for  ammonium  is  unity.  When  ammonium  is  absent,  the 
preference  is  zero.  In  the  presence  of  ammonium  and  nitrate+nitrite,  the 
preference  depends  on  the  abundance  of  both  forms  relative  to  the  half-saturation 
constant  for  nitrogen  uptake.  When  both  ammonium  and  nitrate+nitrite  are 
abundant,  the  preference  for  ammonium  approaches  unity.  When  ammonium  is 
scarce  but  nitrate+nitrite  is  abundant,  the  preference  decreases  in  magnitude  and 
a  significant  fraction  of  algal  nitrogen  requirement  comes  from  nitrate+nitrite. 


Effect  of  algae  on  dissolved  oxygen 

Algae  produce  oxygen  during  photosynthesis  and  consume  oxygen  through 
respiration.  The  quantity  produced  depends  on  the  form  of  nitrogen  utilized  for 
growth.  More  oxygen  is  produced  per  unit  of  carbon  fixed  when  nitrate  is  the 
algal  nitrogen  source  than  when  ammonium  is  the  source.  The  following 
equations  describe  algal  uptake  of  carbon  and  nitrogen  and  production  of 
dissolved  oxygen  (Morel  1983). 
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(3-14) 


1 06  C02  +I6NH4  +H2PC>4  +  106H2O  -» 
protoplasm  +  106O2  +  15H+ 


106CO2  +I6NO3  +H2P04  +  122H20  +  17H+  -> 
protoplasm  + 138  02 


(3-15) 


When  ammonium  is  the  nitrogen  source,  one  mole  oxygen  is  produced  per 
mole  carbon  dioxide  fixed.  When  nitrate  is  the  nitrogen  source,  1.3  moles  oxygen 
are  produced  per  mole  carbon  dioxide  fixed. 

The  following  equation  describes  the  effect  of  algae  on  dissolved  oxygen 
(DO)  in  the  model. 

j{DO  =  [(1. 3-0.3*  PN)*G -(1  - FCD) *  R]*  AOCR *  B  (3-16) 


in  which: 

FCD  =  fraction  of  algal  metabolism  recycled  as  dissolved  organic  carbon 
(0<FCD<  1) 

AOCR  =  dissolved  oxygen-to-carbon  ratio  in  respiration  (2.67  g  02  g'1  C) 

The  magnitude  of  AOCR  is  derived  from  a  simple  representation  of  the 
respiration  process: 

CH20  +  02  =  C02  +  H20  (3-17) 

The  quantity  (1.3  -  0.3  •  PN)  is  the  photosynthesis  ratio  and  expresses  the 
molar  quantity  of  oxygen  produced  per  mole  of  carbon  fixed.  The  photosynthesis 
ratio  approaches  unity  as  the  algal  preference  for  ammonium  approaches  unity. 


Organic  Carbon 

Organic  carbon  undergoes  innumerable  transformations  in  the  water  column. 
The  model  carbon  cycle  (Figure  3-7)  consists  of  the  following  elements: 

•  Phytoplankton  production  and  excretion. 

•  Predation  on  phytoplankton. 

•  Dissolution  of  particulate  carbon. 

•  Heterotrophic  respiration. 

•  Denitrification. 

•  Settling. 
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Figure  3-7.  Model  carbon  cycle 


Algal  production  is  the  primary  carbon  source  although  carbon  also  enters 
the  system  through  external  loading.  Predation  on  algae  by  zooplankton  and 
other  organisms  releases  particulate  and  dissolved  organic  carbon  to  the  water 
column.  A  fraction  of  the  particulate  organic  carbon  undergoes  first-order 
dissolution  to  dissolved  organic  carbon.  Dissolved  organic  carbon  produced  by 
excretion,  by  predation,  and  by  dissolution  is  respired  at  a  first-order  rate  to 
inorganic  carbon.  Particulate  organic  carbon  that  does  not  undergo  dissolution 
settles  to  the  bottom  sediments. 

Dissolution  and  respiration 

Organic  carbon  dissolution  and  respiration  are  represented  as  first-order 
processes  in  w'hich  the  reaction  rate  is  proportional  to  concentration  of  the 
reactant.  An  exponential  function  (Figure  3-5)  relates  dissolution  and  respiration 
to  temperature. 

In  the  model,  aMonod-like  function  (Monod  1949)  diminishes  respiration  as 
dissolved  oxygen  approaches  zero.  As  oxygen  is  depleted  from  natural  systems, 
oxidation  of  organic  matter  is  effected  by  the  reduction  of  alternate  oxidants.  The 
sequence  in  which  alternate  oxidants  are  employed  is  determined  by  the 
thermodynamics  of  oxidation-reduction  reactions.  The  first  substance  reduced  in 
the  absence  of  oxygen  is  nitrate.  A  representation  of  the  denitrification  reaction 
can  be  obtained  by  balancing  standard  half-cell  redox  reactions  (Stumm  and 
Morgan  1981): 

4NO3  +  4H+  +  5CH20  ->  2N2  +  7H20  +  5C02  (3-18) 
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Equation  3-18  describes  the  stoichiometry  of  the  denitrification  reaction.  The 
kinetics  of  the  reaction,  as  represented  in  the  model,  are  first-order.  The  dissolved 
organic  carbon  respiration  rate  Kdoc  is  modified  so  that  significant  decay  via 
denitrification  occurs  only  when  nitrate  is  freely  available  and  dissolved  oxygen 
is  depleted  (Figure  3-8).  A  parameter  is  included  so  that  the  anoxic  respiration 
rate  is  slower  than  oxic  respiration: 


Denit : 

in  which 
Denit  - 
KHodoc -- 

KHndn  - 

AANOX-- 

Kdoc  -- 


KHodoc 


NO3 


KHodoc  +  DO  KHndn  +  NO, 


» AANOX  •  Kdoc 


(3-19) 


denitrification  rate  of  dissolved  organic  carbon  (day'1) 

half-saturation  concentration  of  dissolved  oxygen  required  for  oxic 
respiration  (g  02  m'3) 

half-saturation  concentration  of  nitrate  required  for  denitrification 
(g  N  m'3) 

ratio  of  denitrification  to  oxic  carbon  respiration  rate 
{0  <  AANOX  <  1) 

first-order  dissolved  organic  carbon  respiration  rate  (day4) 


1 

0.9  + 
0.8 
0.7 
0.6 
0.5 
0.4 
0.3  + 
0.2 
0.1 
0 
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N03  /  KHndn 
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Figure  3-8.  Effect  of  dissolved  oxygen  and  nitrate  on  denitrification 


Dissolved  organic  carbon 

The  complete  representation  of  dissolved  organic  carbon  sources  and  sinks  in 
the  model  ecosystem  is 
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—  DOC  =  FCD  •R*B  +  FCDP  •  PR  +  Klpoc  •  LPOC 
St 

+  Krpoc  •  RPOC - — — - •  Kdoc  •  DOC  —  DENIT  •  DOC 

1  KHodoc  +  DO 

in  which 

DOC  =  dissolved  organic  carbon  (g  m'3) 

FCDP  =  fraction  of  predation  on  algae  released  as  DOC  (0  <  FCDP  <  1) 

Klpoc  =  dissolution  rate  of  LPOC  (day-1) 

LPOC  =  labile  particulate  organic  carbon  (g  m'3) 

Krpoc  =  dissolution  rate  of  RPOC  (day'1) 

RPOC  =  refractory  particulate  organic  carbon  (g  m'3) 

Labile  particulate  organic  carbon 

The  complete  representation  of  labile  particulate  organic  carbon  sources  and 
sinks  in  the  model  ecosystem  is 

«  C 

— LPOC  =  FCL  •  R»B  +  FCLP  »PR-  Klpoc  •  LPOC  -Wl»—LPOC  (3-21) 
St  Sz 

in  which 

FCL  =  fraction  of  algal  respiration  released  as  LPOC  (0  <  FCL  <  1) 

FCLP  =  fraction  of  predation  on  algae  released  as  LPOC  (0  <  FCLP  <  1) 

Wl  =  settling  velocity  of  labile  particles  (m  day'1) 

Refractory  particulate  organic  carbon 

The  complete  representation  of  refractory'  particulate  organic  carbon  sources 
and  sinks  in  the  model  ecosystem  is 

— RPOC  =  FCR •R»B  +  FCRP  •  PR-  Krpoc  •  RPOC  - Wr  • — RPOC  (3-22) 
St  Sz 

in  which 

FCR  =  fraction  of  algal  respiration  released  as  RPOC  (0  <  FCR  <  1) 

FCRP  =  fraction  of  predation  on  algae  released  as  RPOC  (0  <  FCRP  <  1) 

Wr  =  settling  velocity  of  refractory  particles  (m  day'1) 
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Phosphorus 

The  model  phosphorus  cycle  (Figure  3-9)  includes  the  following  processes: 

•  Algal  uptake  and  excretion. 

•  Predation. 

•  Hydrolysis  of  particulate  organic  phosphorus. 

•  Mineralization  of  dissolved  organic  phosphorus. 

•  Settling  and  resuspension. 

External  loads  provide  the  ultimate  source  of  phosphorus  to  the  system. 
Dissolved  phosphate  is  incorporated  by  algae  during  growth  and  released  as 
phosphate  and  organic  phosphorus  through  respiration  and  predation.  Dissolved 
organic  phosphorus  is  mineralized  to  phosphate.  A  portion  of  the  particulate 
organic  phosphorus  hydrolyzes  to  dissolved  organic  phosphorus.  The  balance 
settles  to  the  sediments.  Within  the  sediments,  particulate  phosphorus  is 
mineralized  and  recycled  to  the  water  column  as  dissolved  phosphate. 


Figure  3-9.  Model  phosphorus  cycle 


Hydrolysis  and  mineralization 

Within  the  model,  hydrolysis  is  defined  as  the  process  by  which  particulate 
organic  substances  are  converted  to  dissolved  organic  form.  Mineralization  is 
defined  as  the  process  by  w'hich  dissolved  organic  substances  are  converted  to 
dissolved  inorganic  form.  Conversion  of  particulate  organic  phosphorus  to 
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phosphate  proceeds  through  the  sequence  of  hydrolysis  and  mineralization. 

Direct  mineralization  of  particulate  organic  phosphorus  does  not  occur. 

Mineralization  of  organic  phosphorus  is  mediated  by  the  release  of 
nucleotidase  and  phosphatase  enzymes  by  bacteria  (Ammerman  and  Azam  1985; 
Chrost  and  Overbeck  1987)  and  algae  (Matavulj  and  Flint  1987;  Chrost  and 
Overbeck  1987;  Boni  et  al.  1989).  Since  the  algae  themselves  release  the  enzyme 
and  since  bacterial  abundance  is  related  to  algal  biomass,  the  rate  of  organic 
phosphorus  mineralization  is  related,  in  the  model,  to  algal  biomass.  A  most 
remarkable  property  of  the  enzyme  process  is  that  alkaline  phosphatase  activity  is 
inversely  proportional  to  ambient  phosphate  concentration  (Chrost  and  Overbeck 
1987;  Boni  et  al.  1989).  Put  in  different  terms,  when  phosphate  is  scarce,  algae 
stimulate  production  of  an  enzyme  that  mineralizes  organic  phosphorus  to 
phosphate.  This  phenomenon  is  simulated  by  relating  mineralization  to  the  algal 
phosphorus  nutrient  limitation.  Mineralization  is  highest  when  algae  are  strongly 
phosphorus  limited  and  is  least  when  no  limitation  occurs. 

The  expression  for  mineralization  rate  is 

Kdop  =  Kdp  +  — ^ —  •  Kdpalg  •  B  (3-23) 

F  F  KHp  +  P04 

in  which: 

Kdop  =  mineralization  rate  of  dissolved  organic  phosphorus  (day  ') 

Kdp  =  minimum  mineralization  rate  (day1) 

KHp  =  half-saturation  concentration  for  algal  phosphorus  uptake  (g  P  m'3) 
Kdpalg  =  constant  that  relates  mineralization  to  algal  biomass  (m3  g'1  C  day'1) 

Potential  effects  of  algal  biomass  and  nutrient  limitation  on  the 
mineralization  rate  are  shown  in  Figure  3-10.  When  nutrient  concentration 
greatly  exceeds  the  half-saturation  concentration  for  algal  uptake,  the  rate 
roughly  equals  the  minimum.  Algal  biomass  has  little  influence.  As  nutrient 
becomes  scarce  relative  to  the  half-saturation  concentration,  the  rate  increases. 
The  magnitude  of  the  increase  depends  on  algal  biomass.  Factor  of  two  to  three 
increases  are  feasible. 

Exponential  functions  (Figure  3-5)  relate  mineralization  and  hydrolysis  rates 
to  temperature. 
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Figure  3-10.  Effect  of  algal  biomass  and  nutrient  concentration  on  phosphorus 
mineralization 


Dissolved  phosphate 

The  mass-balance  equation  for  dissolved  phosphate  is 
—  PO  4  =  Kdop»DOP-APC»G*B 

St  (3-2- 

+  APC  •[fPI  •  R»  B  +  FPIP  •  PR.] 

in  which 

FPI  =  fraction  of  algal  metabolism  released  as  dissolved  phosphate 
(0  <FPI<  1) 

FPIP  =  fraction  of  predation  released  as  dissolved  phosphate  (0  <FPIP  <  1) 

Dissolved  organic  phosphorus 

The  mass  balance  equation  for  dissolved  organic  phosphorus  is 


— DOP  =  APC  •(R»B»  FPD  +  PR  •  FPDP)  +  Klpop  •  LPOP 
St 

+  Krpop  •  RPOP  -  Kdop  •  DOP 


(3-25) 


in  which: 

FPD  —  fraction  of  algal  metabolism  released  as  DOP  (0  <  FPD  <  1) 
FPDP  =  fraction  of  predation  on  algae  released  as  DOP  (0  <  FPDP  <  1) 
Klpop  =  hydrolysis  rate  of  LPOP  (day'1) 

Krpop  =  hydrolysis  rate  of  RPOP  (day'1) 
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Labile  particulate  organic  phosphorus 

The  mass  balance  equation  for  labile  particulate  organic  phosphorus  is 


S_ 

St 


LPOP  =  APC  •(R»B*  FPL  +  PR  •  FPLP) 


C 

-Klpop  •  LPOP  -  Wl  •  — LPOP 
Sz 


(3-26) 


in  which 

FPL  =  fraction  of  algal  metabolism  released  as  LPOP  (0  <  FPL  <  1) 

FPLP  =  fraction  of  predation  on  algae  released  as  LPOP  (0  <  FPLP  <  1) 

Refractory  particulate  organic  phosphorus 

The  mass  balance  equation  for  refractory  particulate  organic  phosphorus  is 


—RPOP  =  APC  •(R»B»  FPR  +  PR  •  FPRP) 
St 


-Krpop  •  RPOP  -  Wr  •  — RPOP 
Sz 


(3-27) 


in  which 

FPR  =  fraction  of  algal  metabolism  released  as  RPOP  (0  <  FPR  <  1) 

FPRP  =  fraction  of  predation  on  algae  released  as  RPOP  (0  <  FPRP  <  1) 

Nitrogen 

The  model  nitrogen  cycle  (Figure  3-1 1)  includes  the  following  processes: 

•  Algal  production  and  metabolism. 

•  Predation. 

•  Hydrolysis  of  particulate  organic  nitrogen. 

•  Mineralization  of  dissolved  organic  nitrogen. 

•  Settling. 

•  Nitrification. 

•  Denitrification. 

External  loads  provide  the  ultimate  source  of  nitrogen  to  the  system. 
Available  nitrogen  is  incorporated  by  algae  during  growth  and  released  as 
ammonium  and  organic  nitrogen  through  respiration  and  predation.  A  portion  of 
the  particulate  organic  nitrogen  hydrolyzes  to  dissolved  organic  nitrogen.  The 
balance  settles  to  the  sediments.  Dissolved  organic  nitrogen  is  mineralized  to 
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ammonium.  In  an  oxygenated  water  column,  a  fraction  of  the  ammonium  is 
subsequently  oxidized  to  nitrate+nitrite  through  the  nitrification  process.  In 
anoxic  water,  nitrate+nitrite  is  lost  to  nitrogen  gas  through  denitrification. 
Particulate  nitrogen  that  settles  to  the  sediments  is  mineralized  and  recycled  to 
the  water  column,  primarily  as  ammonium.  Nitrate+nitrite  moves  in  both 
directions  across  the  sediment-water  interface,  depending  on  relative 
concentrations  in  the  water  column  and  sediment  interstices. 


Figure  3-11.  Model  nitrogen  cycle 


Nitrification 

Nitrification  is  a  process  mediated  by  specialized  groups  of  autotrophic 
bacteria  that  obtain  energy  through  the  oxidation  of  ammonium  to  nitrite  and 
oxidation  of  nitrite  to  nitrate.  A  simplified  expression  for  complete  nitrification 
(Tchobanoglous  and  Schroeder  1987)  is 

NH4  +  2  02  ->  NO3  +  H20  +  2H+  (3-28) 

The  simplified  stoichiometry  indicates  that  two  moles  of  oxygen  are  required 
to  nitrify  one  mole  of  ammonium  into  nitrate.  The  simplified  equation  is  not 
strictly  true,  however.  Cell  synthesis  by  nitrifying  bacteria  is  accomplished  by  the 
fixation  of  carbon  dioxide  so  that  less  than  two  moles  of  oxygen  are  consumed 
per  mole  ammonium  utilized  (Wezemak  and  Gannon  1968). 

The  kinetics  of  complete  nitrification  are  modeled  as  a  function  of  available 
ammonium,  dissolved  oxygen,  and  temperature: 


3-20 


Chapter  3  Water  Quality  Model  Formulation 


f(T)  •  NTm 


(3-29) 


NT  = 


DO  ^  NH4 
KHont  +  DO  KHnnt  +  NH4 


in  which 

iVT  =  nitrification  rate  (g  N  m'3  day’1) 

KHont  =  half-saturation  constant  of  dissolved  oxygen  required  for  nitrification 
(g  02  m'3) 

KHnnt  =  half-saturation  constant  of  NH*  required  for  nitrification  (g  N  m'3) 
NTm  =  maximum  nitrification  rate  at  optimal  temperature  (g  N  m'3  day'1) 

The  kinetics  formulation  (Figure  3-12)  incorporates  the  products  of  two 
Monod-like  functions  (Monod  1949).  The  first  function  diminishes  nitrification 
at  low  dissolved  oxygen  concentration.  The  second  function  expresses  the 
influence  of  ammonium  concentration  on  nitrification.  When  ammonium 
concentration  is  low,  relative  to  KHnnt ,  nitrification  is  proportional  to  ammonium 
concentration.  For  NH4  «  KHnnt ,  the  reaction  is  approximately  first-order.  (The 
first-order  decay  constant  ~  NTm/KHnnt.)  When  ammonium  concentration  is 
large,  relative  to  KHnnt,  nitrification  approaches  a  maximum  rate.  This 
formulation  is  based  on  a  concept  proposed  by  Tuffey  et  al.  (1974).  Nitrifying 
bacteria  adhere  to  benthic  or  suspended  sediments.  When  ammonium  is  scarce, 
vacant  surfaces  suitable  for  nitrifying  bacteria  exist.  As  ammonium  concentration 
increases,  bacterial  biomass  increases,  vacant  surfaces  are  occupied,  and  the  rate 
of  nitrification  increases.  The  bacterial  population  attains  maximum  density  when 
all  surfaces  suitable  for  bacteria  are  occupied.  At  this  point,  nitrification  proceeds 
at  a  maximum  rate  independent  of  additional  increase  in  ammonium 
concentration. 

The  optimal  temperature  for  nitrification  may  be  less  than  peak  temperatures 
that  occur  in  temperate  waters.  To  allow  for  a  decrease  in  nitrification  at 
superoptimal  temperature,  the  effect  of  temperature  on  nitrification  is  modeled  in 
the  Gaussian  form  of  Equation  3-7. 
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Figure  3-12.  Effect  of  dissolved  oxygen  and  ammonium  concentrations  on 
nitrification  rate 


Effect  of  denitrification  on  nitrate 

The  effect  of  denitrification  on  dissolved  organic  carbon  has  been  described. 
Denitrification  removes  nitrate  from  die  system  in  stoichiometric  proportion  to 
carbon  removal: 


— NO?  =  -ANDC  Denit  DOC  (3-30) 

St  3 


in  which  ANDC  is  the  mass  nitrate-nitrogen  reduced  per  mass  dissolved  organic 
carbon  oxidized  (0.933  g  N  g'1  C). 

Nitrogen  mass  balance  equations 

The  mass-balance  equation  for  nitrogen  state  variables  are  written  by 
summing  all  previously  described  sources  and  sinks: 

a.  Ammonium. 


— NH4  =  ANC  •  [(R  •  FNI  -  PN  •  G)  •  B  +  PR  •  FNIP] 


+Kdon  •  DON  -  NT 


(3-31) 


in  which 

FNI  =  fraction  of  algal  metabolism  released  as  NH4  (0  <  FNI  <  1) 
FNIP  =  fraction  of  predation  released  as  NH4  (0  <  FNIP  <  1) 

Kdon  =  dissolved  organic  nitrogen  mineralization  rate  (day4) 
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b.  Nitrate+nitrite. 


—m23=-ANC*(l-PN)*G»B  +  NT 
St 

-ANDC  •  Denit  •  DOC 


(3-32) 


c.  Dissolved  organic  nitrogen. 

— DON  =  ANC  »(R»B»  FND  +  PR  •  FNDP) 
St 

+  Klpon  •  LPON  +  Krpon  •  RPON  -  Kdon  •  DON 


(3-33) 


in  which: 

FND  =  fraction  of  alga 1  metabolism  released  as  DON  (0  <  FND  <  1) 
FNDP  =  fraction  of  predation  on  algae  released  as  DON  (0  <  FNDP  <  1) 
Klpon  =  hydrolysis  rate  of  LPON  (day'1) 

Krpon  =  hydrolysis  rate  of  RPON  (day'1) 

d.  Labile  particulate  organic  nitrogen 

— LPON  =  ANC  *(R*B»  FNL  +  PR  •  FNLP) 

St 

-  Klpon  •  LPON  -  Wl  •  — LPON 
Sz 


in  which: 

FNL  =  fraction  of  algal  metabolism  released  as  LPON  (0  <  FNL  <  1) 
FNLP  =  fraction  of  predation  on  algae  released  as  LPON  (0  <  FNLP  <  1) 

e.  Refractory  particulate  organic  nitrogen 


St 


RPON  =  ANC  •(R»B»  FNR  +  PR  •  FNRP) 


-  Krpon  •  RPON  -Wr» — RPON 
Sz 
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in  which: 

FNR  =  fraction  of  algal  metabolism  released  as  RPON( 0  <  FNR  <  1) 
FNRP  =  fraction  of  predation  on  algae  released  as  RPON  (0  <  FNRP  <  1) 
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Chemical  Oxygen  Demand 

Chemical  oxygen  demand  is  the  concentration  of  reduced  substances  that  are 
oxidized  through  abiotic  reactions.  In  the  model,  chemical  oxygen  demand 
originates  as  methane  released  from  sediments.  Methane  is  quantified  in  units  of 
oxygen  demand  and  modeled  with  the  following  relationship: 

—COD  = - — - •  Kcod  •  COD  (3-36) 

St  KHocod  +  DO 


in  which 

COD  =  chemical  oxygen  demand  concentration  (g  oxygen-equivalents  m'3) 

KHocod  =  half-saturation  concentration  of  dissolved  oxygen  required  for 
exertion  of  chemical  oxygen  demand  (g  O2  m'3) 

Kcod  =  oxidation  rate  of  chemical  oxygen  demand  (day'1) 

An  exponential  function  (Figure  3-5)  describes  the  effect  of  temperature  on 
exertion  of  chemical  oxygen  demand. 


Dissolved  Oxygen 

Sources  and  sinks  of  dissolved  oxygen  in  the  water  column  (Figure  3-13) 
include  the  following: 

•  Algal  photosynthesis 

•  Atmospheric  reaeration 

•  Algal  respiration 

•  Heterotrophic  respiration 

•  Nitrification 

•  Chemical  oxygen  demand 

Reaeration 

The  rate  of  reaeration  is  proportional  to  the  dissolved  oxygen  deficit  in  model 
segments  that  form  the  air-water  interface: 

-DO  =—»(DOs-DO)  (3-37) 

St  A z 


in  which: 

Kr  =  reaeration  coefficient  (m  day'1) 
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Az  =  model  layer  thickness  (m) 

DOs  =  dissolved  oxygen  saturation  concentration  (g  O2  m'3) 


Figure  3-13.  Dissolved  oxygen  sources  and  sinks 


In  free-flowing  streams,  the  reaeration  coefficient  depends  largely  on 
turbulence  generated  by  bottom  shear  stress  (O’Connor  and  Dobbins  1958).  In 
lakes  and  coastal  waters,  however,  wind  effects  may  dominate  the  reaeration 
process  (O’Connor  1983).  For  Lake  Washington,  a  relationship  for  wind-driven 
gas  exchange  (Hartman  and  Hammond  1985)  was  employed: 


Kr  =  Areaer  •  Rv  •  Wms 1,5 


(3-38) 


in  which 

Areaer  =  empirical  constant  (~0.1) 

Rv  =  ratio  of  kinematic  viscosity  of  pure  water  at  20  °C  to  kinematic 
viscosity  of  water  at  specified  temperature  and  salinity 

Wms  =  wind  speed  measured  at  10  m  above  water  surface  (ms'1) 

Hartman  and  Hammond  (1985)  indicate  Areaer  takes  the  value  0.157.  In  the 
present  model,  Areaer  is  treated  as  a  variable  to  allow  for  effects  of  wind 
sheltering,  for  differences  in  height  of  local  w'ind  observations,  and  for  other 
factors. 

An  empirical  function  (Figure  3-14)  that  fits  tabulated  values  of  Rv  is 
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Rv  =  0.54  +  0. 0233  •  T  -  0.0020  •  S 


(3-39) 


in  which  S  is  salinity  (ppt). 

♦  Table,  S^O 

1-3  ■  Table,  S=35  ::  "2:  TBB- •  T  :y- '3’©';f ST 

12-  Equation,  S=0 
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Figure  3-14.  Computed  and  tabluated  values  of  Rv 


:B..« 


Jit  < 

,  'XK--  : 


|lp||gi 


♦  Table,  S^0 
■  Table,  S=35 
Equation,  S=0 
—  Equation,  S~35 


:  ♦ 
m 


i 


Saturation  dissolved  oxygen  concentration  diminishes  as  temperature  and 
salinity  increase.  The  following  empirical  formula  describes  these  effects  (Genet 
etal.  1974): 


DOs  =  14. 5532  -  0. 3821 7  *T+0. 0054258  •  T2 

-  CL  •(!.  665  _10~4 -5.866  _\0~6  »T  +  9.796  _10~*  •T2 ) 

in  which  CL  is  the  chloride  concentration  (=  salinity/1.80655). 

Mass  balance  equation  for  dissolved  oxygen 

— DO  =  AOCR»[(1.3-0.3»PN)»P-(1-FCD)*R]»B 
St 

-AONT  •  NT - — - •  AOCR  •  Kdoc  •  DOC 

KHodoc  +  DO 

DO  Kr 

—  - — - •  Kcod  •  COD  +  —  •  (DOs  -  DO) 

KHocod  +  DO  A z 

in  which 


(3-40) 


(3-41) 
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AOCR  =  oxygen-to-carbon  mass  ratio  in  production  and  respiration  (=  2.67  g 

o2g-'  C) 

AONT  =  oxygen  consumed  per  mass  ammonium  nitrified  (=  4.33  g  02  g  1  N) 


Temperature 

Computation  of  temperature  employs  a  conservation  of  internal  energy 
equation  that  is  analogous  to  the  conservation  of  mass  equation.  For  practical 
purposes,  the  internal  energy  equation  can  be  written  as  a  conservation  of 
temperature  equation.  The  only  source  or  sink  of  temperature  considered  is 
exchange  with  the  atmosphere.  Atmospheric  exchange  is  considered  proportional 
to  the  temperature  difference  between  the  water  surface  and  a  theoretical 
equilibrium  temperature  (Edinger  et  al.  1974): 

— T  = - — - •(Te-T)  (3-42) 

St  p»Cp»  Az 

in  which 

T  =  water  temperature  (°C) 

KT  =  Heat  exchange  coefficient  (watt  m'2  °C’) 
p  =  density  of  water  (1000  kg  m'3) 

Cp  =  specific  heat  of  water  (4200  watt  s  kg'1  °C"’) 

Te  =  equilibrium  temperature  (°C) 


Inorganic  (Fixed)  Solids 

The  only  kinetics  transformation  of  fixed  solids  is  settling: 

-ISS  =  -Wiss»—ISS  3-43 

St  Sz 

in  which: 

ISS  =  fixed  solids  concentration  (g  m'3) 

Wiss  =  solids  settling  velocity  (m  day'1) 


Parameter  Values 

Model  parameter  evaluation  is  a  recursive  process.  Parameters  are  selected 
from  a  range  of  feasible  values,  tested  in  the  model,  and  adjusted  until 
satisfactory  agreement  between  predicted  and  observed  variables  is  obtained. 
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Ideally,  the  range  of  feasible  values  is  determined  by  observation  or  experiment. 
For  some  parameters,  however,  no  observations  are  available.  Then,  the  feasible 
range  is  determined  by  parameter  values  employed  in  similar  models  or  by  the 
judgement  of  the  modeler.  A  review  of  parameter  values  was  included  in 
documentation  of  the  first  application  of  this  model  (Cerco  and  Cole  1994). 
Parameters  from  the  initial  study  were  refined,  where  necessary,  for  the  present 
model.  A  complete  set  of  parameter  values  is  provided  in  Table  3-2. 


Table  3-2 

Parameters  in  Kinetics  Equations 

Symbol 

Definition 

Value 

Units 

AANOX 

ratio  of  anoxic  to  oxic  respiration 

0.5 

0  <AANOX  <_1 

ANC 

nitrogen-to-carbon  ratio  of  algae 

0.15 

g  N  g'1  C 

AOCR 

dissolved  oxygen-to-carbon  ratio  in 
respiration 

2.67 

g  02  g'1  C 

AONT 

mass  dissolved  oxygen  consumed  per  mass 
ammonium  nitrified 

4.33 

g  o2  g’1  N 

APC 

algal  phosphorus-to-carbon  ratio 

0.0165 

gpg'c 

Areaer 

empirical  constant  in  reaeration  equation 

0.08 

BM 

basal  metabolic  rate  of  algae  at  reference 
temperature  Tr 

0.03 

day*1 

BPR 

base  predation  rate 

0.22  (May -Oct.), 
0.045  otherwise 

day*1 

CChl 

carbon-to-chlorophyll  ratio 

50 

g  C  g'1  chi 

FCD 

fraction  of  dissolved  organic  carbon  produced 
by  algal  metabolism 

0.0 

0< FCD < 1 

FCDP 

fraction  of  dissolved  organic  carbon  produced 
by  predation 

0.6 

0< FCDP< 1 

FCL 

fraction  of  labile  particulate  carbon  produced 
by  algal  metabolism 

0.0 

0  < FCL  <  1 

FCLP 

fraction  of  labile  particulate  carbon  produced 
by  predation 

0.12 

0<  FCLP  <  1 

FCR 

fraction  of  refractory  particulate  carbon 
produced  by  algal  respiration 

0.0 

0<  FCR<  1 

FCRP 

fraction  of  refractory  particulate  carbon 
produced  by  predation 

0.28 

0< FCRP < 1 

FND 

fraction  of  dissolved  organic  nitrogen 
produced  by  algal  metabolism 

0.3 

0  <  FND  <  1 

FNDP 

fraction  of  dissolved  organic  nitrogen 
produced  by  predation 

0.35 

0  <  FNDP  <  1 

FNI 

fraction  of  inorganic  nitrogen  produced  by 
algal  metabolism 

0.6 

0<  FNI  <  1 

FNIP 

fraction  of  inorganic  nitrogen  produced  by 
predation 

0.25 

0<  FNIP<  1 

FNL 

fraction  of  labile  particulate  nitrogen  produced 
by  algal  metabolism 

0.075 

0  <  FNL  <  1 

FNLP 

fraction  of  labile  particulate  nitrogen  produced 
by  predation 

0.12 

0  <  FNLP  <  1 

FNR 

fraction  of  refractory  particulate  nitrogen 
produced  by  algal  metabolism 

0.075 

0< FNR <  1 

FNRP 

fraction  of  refractory  particulate  nitrogen 
produced  by  predation 

0.28 

0  <  FNRP  <  1 

j  (Continued)  || 
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jTable  3-2  (continued) 

Symbol 

Definition 

Value  Units  | 

FPD 

fraction  of  dissolved  organic  phosphorus 
produced  by  algal  metabolism 

0.2 

0<  FPD<  1 

FPDP 

fraction  of  dissolved  organic  phosphorus 
produced  by  predation 

0.2 

0< FPDP < 1 

FPI 

fraction  of  dissolved  inorganic  phosphorus 
produced  by  algal  metabolism 

0.4 

0<  FPI  <  1 

FPIP 

fraction  of  dissolved  inorganic  phosphorus 
produced  by  predation 

0.5 

0  <  FPIP  <  1 

FPL 

fraction  of  labile  particulate  phosphorus 
produced  by  algal  metabolism 

0.2 

0<  FPL<  1 

FPLP 

fraction  of  labile  particulate  phosphorus 
produced  by  predation 

0.09 

0<  FPLP  <  1 

FPR 

fraction  of  refractory  particulate  phosphorus 
produced  by  algal  metabolism 

0.2 

0<  FPR  <  1 

FPRP 

fraction  of  refractory  particulate  phosphorus 
produced  by  predation 

0.21 

0< FPRP < 1 

Kcod 

oxidation  rate  of  chemical  oxygen  demand 

0.1 

day'1 

Kdoc 

dissolved  organic  carbon  respiration  rate 

0.0075 

day'1 

Kdon 

dissolved  organic  nitrogen  mineralization  rate 

0.018 

day"1 

Kdp 

minimum  mineralization  rate  of  dissolved 
organic  phosphorus 

0.12 

day"1 

Kdpalg 

constant  that  relates  mineralization  rate  to 
algal  biomass 

0.2 

m3  g'1  C  day'1 

KHn 

half-saturation  concentration  for  nitrogen 
uptake  by  algae 

0.025 

g  N  m*3 

KHndn 

half-saturation  concentration  of  nitrate 
required  for  denitrification 

0.1 

g  N  m"3 

KHnnt 

half-saturation  concentration  of  NH4  required 
for  nitrification 

1.0 

g  N  m'3 

KHocod 

half-saturation  concentration  of  dissolved 
oxygen  required  for  exertion  of  COD 

0.5 

g  02  m'3 

KHodoc 

half-saturation  concentration  of  dissolved 
oxygen  required  for  oxic  respiration 

0.5 

g  02  m'3 

KHont 

half-saturation  concentration  of  dissolved 
oxygen  required  for  nitrification 

3.0 

g  02  m'3 

KHp 

half-saturation  concentration  for  phosphorus 
uptake  by  algae 

0.005 

g  P  m'3 

Klpoc 

labile  particulate  organic  carbon  dissolution 
rate 

0.005 

day'1 

Klpon 

labile  particulate  organic  nitrogen  hydrolysis 
rate 

0.08 

day'1 

Klpop 

labile  particulate  organic  phosphorus 
hydrolysis  rate 

0.1 

day'1 

Krpoc 

refractory  particulate  organic  carbon 
dissolution  rate 

0.001 

day'1 

Krpon 

refractory  particulate  organic  nitrogen 
hydrolysis  rate 

0.001 

day'1 

Krpop 

refractory  particulate  organic  phosphorus 
hydrolysis  rate 

0.001 

day'1 

KTb 

effect  of  temperature  on  basal  metabolism  of 
algae 

0.032 

"C'1 

j  (Continued)  || 
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Table  3-2  (concluded) 

Symbol 

Definition 

Value 

Units 

KTcod 

effect  of  temperature  on  exertion  of  chemical 
oxygen  demand 

0.041 

day1 

KTgl 

effect  of  temperature  below  Topt  on  growth  of 
algae 

0.003 

°c-2 

KTg2 

effect  of  temperature  above  Topt  on  growth  of 
algae 

0.01 

°c-2 

KThdr 

effect  of  temperature  on  hydrolysis  rates 

0.069 

°C'1 

KTmnl 

effect  of  temperature  on  mineralization  rates 

0.069 

°c-i 

KTntl 

effect  of  temperature  below  Tmnt  on 
nitrification 

0.001 

°C-2 

KTnt2 

effect  of  temperature  above  Tmnt  on 
nitrification 

0.001 

Oq-2 

NTm 

maximum  nitrification  rate  at  optimal 
temperature 

0.075 

g  N  m'3  day'1 

PBm 

maximum  photosynthetic  rate 

250 

g  C  g'1  Chi  day'1 

Presp 

photorespiration  fraction 

0.25 

0<  Presp  <1 

Tmnt 

optimal  temperature  for  nitrification 

30 

°C 

Topt 

optimal  temperature  for  growth  of  algae 

25 

°C 

Tref 

reference  temperature  for  metabolism 

20 

°C 

Trhdr 

reference  temperature  for  hydrolysis 

20 

°c 

Trmnl 

reference  temperature  for  mineralization 

20 

°c 

Wa 

algal  settling  rate 

0.1 

m  day"1 

Wl 

settling  velocity  of  labile  particles 

0.8 

m  day*1 

Wiss 

settling  velocity  of  fixed  solids 

1.0 

m  day'1 

Wr 

settling  velocity  of  refractory  particles 

0.8 

m  day"1 

a 

initial  slope  of  production  versus  irradiance 
relationship 

8.0 

GCcf1  Chi 
(E  n-ry 
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4  Pathogens  and  Toxicants 


Introduction 

At  the  request  of  the  sponsor,  a  pathogen  and  two  toxicants  were  added  to  the 
model  code.  The  pathogen  is  activated  in  the  present  application  and  represents 
fecal  coliform  bacteria.  The  two  toxicants  are  coded  in  a  generalized  format  and 
are  not  presently  activated.  The  general  formulation  is  intended  to  provide  a  rapid 
screening  tool.  The  formulation  also  provides  the  foundation  for  a  more  detailed 
application,  should  this  be  desired. 


Conservation  of  Mass  Equation 

The  foundation  of  CE-QUAL-ICM  is  the  solution  to  the  three-dimensional 
mass-conservation  equation  for  a  control  volume.  Control  volumes  correspond  to 
cells  on  the  model  grid.  CE-QUAL-ICM  solves,  for  each  volume  and  for  each 
state  variable,  the  equation: 


8  Vj  •  Cj 


8  t 


=  TQk  •  c*  +  £a  •  Dk  • 

k=l  k=l 


8  C 


XJc 


+  2  Sj 


(4-1) 


in  which: 


Vj  -  volume  of fh  control  volume  (m3) 

Cj  =  concentration  in  fh  control  volume  (g  m'3) 

t,  x  =  temporal  and  spatial  coordinates 

n  =  number  of  flow  faces  attached  to  fh  control  volume 

Qk  =  volumetric  flow  across  flow  face  k  of fh  control  volume  (m3  s'1) 

Ck  =  concentration  in  flow  across  face  k  (g  m'3) 

Ak=  area  of  flow  face  k  (m2) 

Dk  =  diffusion  coefficient  at  flow  face  k  (m2  s'1) 

Sj  =  external  loads  and  kinetic  sources  and  sinks  in  fh  control  volume 
(gs1) 
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Solution  of  Equation  4-1  on  a  digital  computer  requires  discretization  of  the 
continuous  derivatives  and  specification  of  parameter  values.  The  equation  is 
solved  using  the  QUICKEST  algorithm  (Leonard  1979)  in  the  horizontal  plane 
and  a  Crank-Nicolson  scheme  in  the  vertical  direction.  Discrete  time-steps, 
determined  by  computational  stability  requirements,  are  approximately 
15  minutes. 

The  remainder  of  this  chapter  is  devoted  to  detailing  the  kinetics  sources  and 
sinks  and  to  reporting  parameter  values.  For  notational  simplicity,  the  transport 
terms  are  dropped  in  the  reporting  of  kinetics  formulations. 


Pathogen 

The  pathogen  undergoes  temperature -dependent  first-order  decay: 


S  t 


PATH  =  -  Kpath  •  f(T)  •  PATH 


(4-2) 


in  which: 

PATH  =  pathogen  concentration  (mpn  100  mL'1) 
Kpath  =  decay  rate  at  reference  temperature  T  (day-1) 


For  Lake  Washington,  the  pathogen  decay  rate  is  0. 1  day'1  at  20  °C.  The 
temperature  function,/^!),  is  an  exponential  relationship  in  which  decay  rate 
doubles  for  a  10°  increase  in  temperature  (Figure  4-1). 


Chapter  4  Pathogens  and  Toxicants 


Toxicants 


The  two  toxicants  are  subject  to  identical  kinetics  processes  (Figure  4-2) 
including  the  following: 

•  Decay. 

•  Volatilization. 

•  Partitioning  to  solids. 

•  Settling. 

•  Burial. 


k  volatilization  (Kvol) 


air-water 

interface 


(Ka  dtox) 

Fd  -◄ - ►  Fp 

exchange  with  solids 


decay  (Ktox) 


settling  of 

diffusion  of 

particulate  fraction 

^  dissolved  fraction 

sediment-water  (Wnet) 

(s) 

interface 

(Kadsedtox) 

if  Fd  < - ►Fp 

decay  (Ktoxsed) 


exchange  with  solids 


^  burial  to  deep  sediments  (Wbur) 


Figure  4-2.  Schematic  of  toxicant  model 


Toxicant  1  partitions  to  inorganic  solids  (which  can  be  configured  to 
represent  a  metal).  Toxicant  2  partitions  to  particulate  organic  carbon.  Toxicant 
not  lost  through  decay  or  volatilization  may  settle  to  the  bottom  sediments. 
Within  the  sediments,  toxicant  is  subject  to  decay  and  burial.  Dissolved  toxicant 
diffuses  in  either  direction  across  the  sediment-water  interface. 


Reactions  in  the  water  column 

The  basic  representations  for  the  two  toxicants  are  identical: 

c  TOX 

—  TOX  =  - Ktox  •  f(T)  •  TOX -Kvol  •  Fd • - 

S  t  A z 

-Ws  •  jr  (Fp»TOX) 


in  which: 
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TOX  =  toxicant  concentration  (g  m'3) 

Ktox  =  decay  rate  at  reference  temperature  T  (day'1) 

Kvol  =  volatilization  rate  (m  day'1) 

Fd  =  dissolved  fraction  of  total  toxicant  (0  <  Fd  <  1) 

Fp  =  particulate  fraction  of  total  toxicant  (=  1  -  Fd) 
fVs  =  particle  settling  velocity  (m  day'1) 

A z  =  surface  layer  thickness  (m) 

The  decay  rate  increases  as  an  exponential  function  of  temperature.  The 
volatilization  formulation  assumes  that  toxicant  concentration  in  the  atmosphere 
is  negligible. 

Particulate  and  dissolved  fractions 

For  Toxicant  1,  the  particulate  fraction  is: 


KADtox  1 » ISS 
1  +  KADtox]  •  ISS 


(4-4) 


in  which: 

KADtox  1  =Toxicant  1  partition  coefficient  (m3  g'1) 

ISS  =  inorganic  solids  concentration  (g  m'3) 

For  Toxicant  2,  the  particulate  fraction  is: 

KADtox 2  *  (B  +  LPOC  +  RPOC) 

P  ~  1  +  KADtox 2  •  (B  +  LPOC  +  RPOC ) 

in  which: 

KADtox2  =Toxicant  2  partition  coefficient  (m3  g'1  C) 

B  =  algal  biomass  (g  C  m'3) 

LPOC  =  labile  particulate  organic  carbon  (g  C  m'3) 
RPOC  =  refractory  particulate  organic  carbon  (g  C  m'3) 

The  dissolved  fraction  for  both  toxicants  is: 


Fd  =  \-Fp 


(4-6) 
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Settling 


The  appropriate  settling  velocity  for  Toxicant  1  is  the  inorganic  solids 
settling  velocity,  Wiss.  The  settling  velocities  for  the  multiple  components  to 
which  Toxicant  2  partitions  may  vary.  The  settling  term  is  weighted  by  the 
product  of  the  component  settling  velocities  and  concentrations: 


Ws»~(Fp*TOX)  = 

8  ts- a  r\  .  WSa  »B  +  WSl •  LPOC  +  WSr  •  RPOC 

— KADtoxl  • - 7 - — — r-  •  1  uX 

Sz  1  +  KADtoxl  »(B  +  LPOC  +  RPOC) 


in  which: 

WSa  =  algal  settling  velocity  (m  day'1) 

WSl  =  labile  particulate  organic  carbon  settling  velocity  (m  day'1) 

WSr  =  refractory  particulate  organic  carbon  settling  velocity  (m  day'1) 

Reactions  in  the  sediments 

Bed  sediments  are  envisioned  as  a  single,  well-mixed  layer  (Figure  4-2). 
Toxicants  are  exchanged  with  the  overlying  water  through  settling  of  the 
particulate  fraction  and  diffusion  of  the  dissolved  fraction.  Within  the  sediments, 
toxicants  undergo  decay  and  burial  to  deep,  inactive  sediments.  The  mass- 
balance  equations  for  both  toxicants  are  identical: 

H  •  —TOXsed  =  -( Ktoxsed  •  f  (T)  •  H  +  Wbur)  •  TOXsed 

St  (4-8 

+  Wnet  •  Fpw  •  TOXw  -  s  •  ( Fdsed  •  TOXsed  -  Fdw  •  TOXw ) 

in  which: 

H= 

TOXsed  = 

Ktoxsed = 

Wbur  = 

Wnet  = 

Fpw  = 

TOXw  = 
s  = 

Fdsed  = 

Fdw  = 


thickness  of  active  sediment  layer  (m) 

bulk  concentration  of  toxicant  in  sediments  (g  m'3) 

decay  rate  in  sediments  (day'1) 

burial  rate  to  deep,  inactive  sediments  (m  day'1) 

net  settling  velocity  of  solids  into  sediments  (m  day'1) 

particulate  fraction  in  water  (=  1  -  Fdw ) 

concentration  of  toxicant  in  overlying  water  (g  m'3) 

sediment-water  mass-transfer  velocity  (m  day'1) 

dissolved  fraction  in  sediments  (0  <  Fdsed  <  1) 

dissolved  fraction  in  the  overlying  water  (0  <  Fdw  <  1) 
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Decay  rates  and  partition  coefficients  within  the  sediments  may  vary  from 
those  specified  for  the  water  column.  For  Toxicant  1,  the  dissolved  fraction 
within  the  sediments  is 


Fdsed  = 


_ 1 _ 

1  +  KADsedtox  1  •  ISSsed 


(4-9) 


in  which: 

KADsedtoxl  =  Toxicant  1  partition  coefficient  in  sediments  (m3  g'1) 
ISSsed  =  bulk  inorganic  solids  concentration  in  sediments  (g  rrf3) 

For  Toxicant  2,  the  dissolved  fraction  is 


Fdsed  = 


_ 1 _ 

1  +  KADsedtoxl  •  (G1  +  G2  +  G3) 


(4-10) 


In  which 

KADsedtox2  =  Toxicant  2  partition  coefficient  in  sediments  (m3  g'1  C) 

G1  =  labile  organic  carbon  concentration  in  sediments  (g  C  m'3) 

G2  =  refractory  organic  carbon  concentration  in  sediments  (g  C  m'3) 
G3  =  inert  organic  carbon  concentration  in  sediments  (g  C  m'3) 

The  bulk  sediment  solids  concentration  is  an  input  to  the  diagenesis  model. 
The  sediment  carbon  concentrations  are  computed  within  the  sediment  model. 
The  sediment-water  mass-transfer  velocity  is  an  input  to  the  toxicant  model.  A 
small  amount  of  recoding  would  allow  the  use  of  the  mass-transfer  velocity 
computed  in  the  diagenesis  model. 


Reference 

Leonard,  B.  (1979).  “A  stable  and  accurate  convection  modelling  procedure 
based  on  quadratic  upstream  interpolation,”  Computer  Methods  in  Applied 
Mechanics  and  Engineering  19,  59-98. 
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5  Modeling  Processes  at  the 
Sediment-Water  Interface 


Introduction 

The  predictive  benthic  sediment  model  applied  to  Lake  Washington  was  first 
developed  for  use  in  Chesapeake  Bay  (DiToro  and  Fitzpatrick  1993). 
Management  of  the  bay  required  a  model  with  two  fundamental  capabilities: 

•  Predict  effects  of  management  actions  on  sediment-water  exchange 
processes. 

•  Predict  time  scale  for  alterations  in  sediment-water  exchange  processes. 

The  model  (Figure  5-1,  Table  5-1)  is  driven  by  net  settling  of  organic  matter 
from  the  water  column  to  the  sediments.  In  the  sediments,  the  model  simulates 
the  diagenesis  (decay)  of  the  organic  matter.  Diagenesis  produces  oxygen 
demand  and  inorganic  nutrients.  Oxygen  demand,  as  sulfide  (in  salt  water)  or 
methane  (in  fresh  water),  takes  three  paths  out  of  the  sediments:  export  to  the 
water  column  as  chemical  oxygen  demand,  oxidation  at  the  sediment-water 
interface  as  sediment  oxygen  demand,  or  burial  to  deep,  inactive  sediments. 
Inorganic  nutrients  produced  by  diagenesis  take  two  paths  out  of  the  sediments: 
release  to  the  water  column  or  burial  to  deep,  inactive  sediments. 

The  formulation  of  the  diagenesis  model  is  too  extensive  to  repeat  here. 
Complete  model  documentation  was  provided  by  DiToro  and  Fitzpatrick  (1993). 
More  accessible  documentation  may  presently  be  found  in  DiToro  (2001). 

Details  of  the  sediment  model  required  to  understand  the  coupling  of  the 
sediment  submodel  to  the  model  of  the  water  column  are  provided  in  this  chapter. 
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Figure  5-1 .  Sediment  model  schematic 


Table  5-1 

Sediment  Model  State  Variables  and  Fluxes 

State  Variable 

Sediment-Water  Flux 

Temperature 

Particulate  Organic  Carbon 

Sediment  Oxygen  Demand 

Sulfide/Methane 

Release  of  Chemical  Oxygen  Demand 

Particulate  Organic  Nitrogen 

Ammonium 

Ammonium  Flux 

Nitrate 

Nitrate  Flux 

Particulate  Organic  Phosphorus 

Phosphate 

Phosphate  Flux 

Particulate  Biogenic  Silica1 

Available  Silica1 

Silica  Flux1 

| 1  Not  activated  in  Lake  Washington.  | 

Coupling  With  the  Sediment  Diagenesis  Model 

Benthic  sediments  are  represented  as  two  layers  with  a  total  depth  of  10  cm 
(Figure  5-2).  The  upper  layer,  in  contact  with  the  water  column,  may  be  oxic  or 
anoxic  depending  on  dissolved  oxygen  concentration  in  the  water.  The  lower 
layer  is  permanently  anoxic.  The  thickness  of  the  upper  layer  is  determined  by 
the  penetration  of  oxygen  into  the  sediments.  At  its  maximum  thickness,  the  oxic 
layer  depth  is  only  a  small  fraction  of  the  total. 
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Figure  5-2.  Sediment  model  elevation 


The  sediment  model  consists  of  three  basic  processes.  The  first  is  deposition 
of  particulate  organic  matter  from  the  water  column  to  the  sediments.  Due  to  the 
negligible  thickness  of  the  upper  layer,  deposition  proceeds  from  the  water 
column  directly  to  the  lower,  anoxic  layer.  Within  the  lower  layer,  organic  matter 
is  subject  to  the  second  basic  process,  diagenesis  (or  decay).  The  third  basic 
process  is  flux  of  substances  produced  by  diagenesis  to  the  upper  sediment  layer, 
to  the  water  column,  and  to  deep,  inactive  sediments.  The  flux  portion  of  the 
model  is  the  most  complex.  Computation  of  flux  requires  consideration  of  the 
following: 

•  Reactions  in  both  sediment  layers. 
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•  Partitioning  between  particulate  and  dissolved  fractions  in  both  layers. 

•  Sedimentation  from  the  upper  to  lower  layer  and  from  the  lower  layer  to 
deep  inactive  sediments. 

•  Particle  mixing  between  layers. 

•  Diffusion  between  layers. 

•  Mass  transfer  between  the  upper  layer  and  the  water  column. 

The  water  quality  and  sediment  models  interact  on  a  time  scale  equal  to  the 
integration  time-step  of  the  water  quality  model.  After  each  integration,  predicted 
particle  deposition,  temperature,  nutrient,  and  dissolved  oxygen  concentrations 
are  passed  from  the  water  quality  model  to  the  sediment  model.  The  sediment 
model  computes  sediment-water  fluxes  of  dissolved  nutrients  and  oxygen  based 
on  predicted  diagenesis  and  concentrations  in  the  sediments  and  water.  The 
computed  sediment-water  fluxes  are  incorporated  by  the  water  quality  model  into 
appropriate  mass  balances  and  kinetic  reactions. 


Deposition 

Deposition  is  one  process  that  couples  the  model  of  the  water  column  with 
the  model  of  the  sediments.  Consequently,  deposition  is  represented  in  both  the 
sediment  and  water-column  models.  In  the  water  column,  deposition  is 
represented  with  a  modification  of  the  mass-balance  equation  applied  only  to 
cells  that  interface  the  sediments: 

—  =  [transport]  +  [kinetics]  +  •  Cup  ~  •  C  (5-1) 

St  Az  Az 


in  which: 

C  =  concentration  of  particulate  constituent  in  cell  above  sediments 
(gm-3) 

WS  =  settling  velocity  in  water  column  (m  day'1) 

Az  =  cell  thickness  (m) 

Cup  =  constituent  concentration  two  cells  above  sediments  (g  m'3) 

Wmt=  net  settling  to  sediments  (m  day'1) 

Net  settling  to  the  sediments  may  be  less  than  or  equal  to  settling  in  the  water 
column.  Sediment  resuspension  is  implied  when  settling  to  the  sediments  is  less 
than  settling  through  the  water  column. 


Diagenesis 

Organic  matter  in  the  sediments  is  divided  into  three  G  classes  or  fractions,  in 
accordance  with  principles  established  by  Westrich  and  Bemer  (1984).  Division 
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into  G  classes  accounts  for  differential  decay  rates  of  organic  matter  fractions. 

The  G1  or  labile  fraction  has  a  half-life  of  20  days.  The  G2  or  refractory  fraction 
has  a  half-life  of  1  year.  The  G3  or  inert  fraction  undergoes  no  significant  decay 
before  burial  into  deep,  inactive  sediments.  Each  G  class  has  its  own  mass- 
conservation  equation: 

H  •  —  =  Wr*t  •  fi  •  C-W  •  Gi-H  •  Ki  •  Gi  •  tf'X)  (5-2) 

St 


in  which: 


H  =  total  thickness  of  sediment  layer  (m) 

Gi  =  concentration  of  organic  matter  in  G  class  i  (g  m'3) 
fi  =  fraction  of  deposited  organic  matter  assigned  to  G  class  i 
W  =  burial  rate  (m  day-1) 

Ki  =  decay  rate  of  G  class  /  (day1) 

0;  =  constant  that  expresses  effect  of  temperature  on  decay  of  G  class  i 


Since  the  G3  class  is  inert,  K3  =  0. 


Sediment-water  flux 

The  exchange  of  dissolved  substances  between  the  sediments  and  water 
column  is  driven  by  the  concentration  difference  between  the  surface  sediment 
layer  and  the  overlying  water.  Flux  may  be  in  either  direction  across  the 
sediment-water  interface,  depending  on  concentration  gradient.  Sediment-water 
flux  is  computed  within  the  diagenesis  model  as  the  product  of  concentration 
difference  and  an  internally  computed  mass-transfer  coefficient.  In  the  water 
column,  sediment-water  exchange  of  dissolved  substances  is  represented  with  a 
modification  of  the  mass-balance  equation  applied  only  to  cells  that  interface 
with  bottom  sediments: 

SC  7  0-7  BENFLX 

—  =  [transport]  +  [kinetics]  + -  \?-3) 

St  A z 

in  which  BENFLX  is  the  sediment-water  flux  of  dissolved  substance 
(g  m‘2  day’1). 

By  convention,  positive  fluxes  are  from  sediment  to  water.  Negative  fluxes, 
including  sediment  oxygen  demand,  are  from  water  to  sediments. 
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Parameter  Specification 

Coupling  with  the  sediment  model  requires  specification  of  net  settling  rates, 
of  the  G  splits  of  organic  matter,  and  of  burial  rates. 


Net  settling  rates 

Net  settling  for  inorganic  solids  and  detritus  was  specified  as  half  the  settling 
rate  in  the  water  column.  The  lower  net  settling  rates,  0.5  m  day'1  for  solids  and 
0.4  m  day'1  for  detritus,  were  specified  in  an  attempt  to  reproduce  high  particle 
concentrations  observed  near  the  bottom  of  the  water  column.  Net  settling  for 
algae  was  the  same  as  settling  through  the  water  column,  0.1  m  day'1. 

Assignment  to  G  classes 

Upon  deposition  of  particulate  organic  matter  in  the  sediments,  state 
variables  representing  particulate  organic  matter  in  the  water  quality  model 
required  conversion  into  sediment  model  state  variables.  The  water  quality  model 
considered  two  classes  of  particulate  organic  matter:  labile  and  refractory.  The 
sediment  model  was  based  on  three  classes  of  organic  particles:  labile  (Gl), 
refractory  (G2),  and  inert  (G3).  Labile  particles  from  the  water  quality  model 
were  transferred  directly  into  the  Gl  class  in  the  sediment  model.  Refractory 
particles  from  the  water  quality  model  had  to  be  split  into  G2  and  G3  fractions 
upon  entering  the  sediments.  Algae  settling  directly  to  the  sediments  also 
required  routing  into  sediment  model  state  variables.  Guidance  for  the  splits  was 
obtained  from  phytoplankton  decomposition  experiments  (Westrich  and  Berner 
1984).  Planktonic  particulate  organic  carbon  was  found  to  be  50  percent  labile, 

16  percent  refractory,  and  34  percent  non-reactive  (Table  5-2). 


Table  5-2 

Routing  Organic  Particles  into  Sediment  Classes 

WQM  Variable 

Carbon 

Nitrogen 

Phosphorus  | 

%  Gl 

%  G2 

%  G3 

%  Gl 

%  G2 

%  G3 

%  Gl 

%  G2 

%  G3 

Labile  Particles 

100 

100 

100 

Refractory  Particles 

32 

68 

32 

68 

32 

68 

Algae 

50 

16 

34 

50 

16 

34 

50 

16 

34 

Burial  rates 

The  burial  rate  in  the  sediment  model,  as  applied  to  central  Chesapeake  Bay, 
is  0.25  cm  year'1.  Kuivila  and  Murray  (1984)  cite  sedimentation  rates  of  0.31  to 
0.53  cm  year'1  in  Lake  Washington.  In  view  of  the  similarity  in  burial  rates  and 
the  minor  role  of  this  rate  in  the  model,  the  burial  rate  of  0.25  cm  year'1  was 
retained  for  this  application. 
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Sediment  model  parameters 

With  two  exceptions,  all  parameters  within  the  diagenesis  model  are  exactly 
as  derived  for  the  original  model  application  (DiToro  and  Fitzpatrick  1993).  The 
exceptions  are  the  phosphate  partition  coefficients  in  the  sufficial  and  subsurface 
sediments.  The  original  values  were  halved  for  the  Lake  Washington  application. 
Phosphorus  partition  coefficients  were  reduced  from  30,000  L  kg'1  to  7,500  L 
kg'1  in  the  surface  sediments  and  from  100  L  kg'1  to  50  L  kg'1  in  the  subsurface 
sediments.  In  applying  the  diagenesis  model  to  multiple  systems,  DiToro  (2001) 
noted  the  phosphorus  partition  coefficients  frequently  required  revision.  The 
authors  have  had  the  same  experience  with  applications  in  the  present  study.  The 
partition  coefficients  employed  in  Lake  Washington  are  low  compared  with  the 
range  reported  by  DiToro.  The  low  partition  coefficients  imply  Lake  Washington 
sediments  have  low  phosphorus  retention  relative  to  the  other  systems.  The 
reduction  of  partition  coefficients  is  the  opposite  of  the  authors’  experience  to 
date.  The  partition  coefficients  for  freshwater  systems  are  usually  increased 
relative  to  the  original  estuarine  application.  The  authors  cannot  speculate  on  the 
apparent  difference  in  properties  for  Lake  Washington. 
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6  Coupling  with  the 
Hydrodynamic  Model 


Introduction 

Modeling  the  physics,  chemistry,  and  biology  of  Lake  Washington  required  a 
pair  of  models.  Transport  processes  were  modeled  by  a  three-dimensional 
hydrodynamic  model  that  operated  independently  of  the  water  quality  model. 
Transport  information  from  the  hydrodynamic  model  was  processed  and  stored 
on-line  for  subsequent  use  by  the  water  quality  model. 


The  Hydrodynamic  Model 

CH3D-WES 

The  CH3D-WES  (Computational  Hydrodynamics  in  Three  Dimensions  - 
Waterways  Experiment  Station)  hydrodynamic  model  was  a  substantially  revised 
version  of  the  CH3D  model  originally  developed  by  Sheng  (1986).  Model 
formulation  was  based  on  principles  expressed  by  the  equations  of  motion, 
conservation  of  volume,  and  conservation  of  mass.  Quantities  computed  by  the 
model  included  three-dimensional  velocities,  surface  elevation,  vertical  viscosity 
and  diffusivity,  temperature,  and  density. 

Computational  grid 

The  basic  equations  of  CH3D-WES  were  solved  via  the  finite-difference 
method.  The  finite-difference  solution  algorithm  replaced  continuous  derivatives 
in  the  governing  differential  equations  with  ratios  of  discrete  quantities.  Solutions 
to  the  hydrodynamics  were  obtained  using  1 -minute  intervals  for  the  discrete 
time-steps.  The  spatial  continuum  of  the  lake  was  divided  into  a  grid  of  discrete 
cells.  To  achieve  close  conformance  of  the  grid  to  lake  geometry,  cells  were 
represented  in  curvilinear  rather  than  rectangular  coordinates.  A  z-plane  grid  was 
employed  in  which  the  number  of  vertical  layers  varied  depending  on  local 
depth.  Velocity  and  diffusivity  were  computed  on  the  boundaries  between  cells. 
Temperature  and  density  were  computed  at  the  center  of  each  cell. 

The  computational  grid  extended  from  the  mouth  of  the  Cedar  River  to  the 
mouth  of  the  Sammamish  River,  along  the  north-south  axis,  and  westward  to  the 
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Chittenden  Locks  (Figure  6-1).  The  grid  contained  655  cells,  roughly  0.15  km2  in 
area,  in  the  surface  plane.  Total  number  of  cells  in  the  grid  was  12,177.  Number 
of  cells  in  the  vertical  ranged  from  1  to  41  (Figures  6-2,  6-3).  Surface  cells  were 
2.14  m  thick  at  mean  surface  level.  Variations  in  surface  level  caused  by  wind, 
runoff,  and  other  factors  were  represented  by  computed  variations  in  thickness  of 
the  surface  layer.  Thickness  of  all  subsurface  cells  was  fixed  at  1.53  m. 

Calibration  and  verification 

The  hydrodynamic  model  was  calibrated  and  verified  primarily  through 
comparison  of  computed  and  observed  temperatures.  Qualitative  examinations  of 
circulation  and  other  factors  were  also  conducted.  Details  of  the  calibration 
procedure  were  presented  by  Johnson  et  al.  (2003). 


Linkage  to  the  Water  Quality  Model 

Hydrodynamics  for  employment  in  the  water  quality  model  were  produced 
for  three  years,  1995  -  1997.  Each  year  was  a  single,  continuous  production  run. 
Initial  temperature  for  the  first  year  was  derived  from  observations.  The  initial 
hydrodynamic  field  was  obtained  from  a  5-day  spin-up  period.  Thereafter,  initial 
conditions  for  each  year  were  taken  from  conditions  computed  at  the  end  of  the 
previous  year.  Consequently,  the  hydrodynamic  simulation  was  effectively  a 
continuous  3-year  run,  initialized  only  once. 

Computed  flows,  surface  elevations,  and  vertical  diffusivities  were  output  at 
1-hour  intervals  for  use  in  the  water  quality  model.  The  1-hour  hydrodynamics 
were  determined  as  arithmetic  means  of  hydrodynamics  computed  on  a  1 -minute 
basis.  The  algorithms  and  codes  for  linking  the  hydrodynamic  and  water  quality 
models  were  developed  over  a  decade  ago  and  have  been  tested  and  proved  in 
extensive  applications  since  then.  In  every  model  application,  the  linkage  is 
verified  by  comparing  transport  of  a  conservative  substance  computed  in  each 
model.  For  this  application,  temperatures  computed  by  the  two  models  were 
compared.  In  each  model,  daily  average  temperature  was  sampled  at  locations 
corresponding  to  sample  stations  in  Lakes  Washington  and  Union.  Eighteen 
stations  were  sampled  once  or  twice  a  month,  from  surface  to  bottom,  over  the  3- 
year  simulation  period.  A  total  of  3,468  comparisons  were  made  between  the  two 
models  and  summarized  in  the  form  of  a  scatterplot  (Figure  6-4). 

Correspondence  between  the  two  models  is  excellent.  The  linear  correlation 
coefficient  of  the  independently  computed  temperatures  is  -0.99  and  relative 
error  is  less  than  1  percent.  Temperatures  computed  by  the  two  models  are 
typically  within  one  degree  of  each  other.  A  few  larger  differences  exist.  The 
authors  attribute  these  primarily  to  differences  in  the  specification  of  runoff 
temperatures.  Minor  discrepancies  may  also  be  attributed  to  specification  of 
initial  conditions  and  to  methods  used  to  obtain  daily  average  temperature.  No 
doubt  exists  that  transport  processes  in  the  water  quality  model  are  identical  to 
those  in  the  hydrodynamic  model. 
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Figure  6-1 .  Plan  view  of  Lake  Washington  grid 
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Figure  6-3.  Plan  view  showing  number  of  vertical  cells 
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Figure  6-4.  Comparison  of  daily  average  temperatures  computed  by 
hydrodynamic  (CH3D)  and  water  quality  (ICM)  models 
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7  Observation-based 
Calibration 


The  Database 

The  primary  database  for  calibration  of  the  model  was  derived  from  a 
monitoring  program  conducted  by  King  County  Department  of  Parks  and  Natural 
Resources  (KCDPNR).  Observations  collected  at  thirteen  stations  in  Lake 
Washington  (Figure  7-1)  and  five  stations  in  Lake  Union  (Figure  7-2)  were 
provided  in  May  2001.  Sampling  surveys  were  conducted  once  or  twice  per 
month  and  sometimes  extended  over  a  multiday  period.  For  most  constituents, 
observations  were  collected  at  multiple  depths  in  the  water  column.  Depth 
interval  varied.  Usually  a  sample  was  collected  at  the  1-m  depth.  Thereafter, 
samples  were  collected  at  intervals  of  4  to  20  m  with  smaller  increments  closer  to 
the  surface.  Fecal  coliform  were  sampled  at  the  surface  only.  For  chlorophyll,  the 
predominant  observations  consisted  of  a  single  sample  integrated  over  the  Secchi 
depth.  At  one  deep  site.  Station  852,  vertical  chlorophyll  profiles  were 
occasionally  collected. 

Observed  quantities  relevant  to  the  model  are  listed  in  Table  7-1.  Total 
chlorophyll  and  phaeophytin  observations  were  reported.  Phaeophytin  was 
subtracted  from  total  chlorophyll  to  obtain  viable  chlorophyll  for  comparison 
with  the  modeled  quantity.  Observed  disk  visibility  was  related  to  modeled  light 
attenuation  through  a  relationship  attributed  to  Poole  and  Atkins  (Walker  1982): 

Ke  ■  DV  =  k  (7-1) 


in  which: 

Ke  =  diffuse  attenuation  coefficient  (m'1) 

DV=  disk  visibility  (m) 
k  =  empirical  constant 

The  value  1.86  was  assigned  to  k  based  on  the  median  value  observed  in 
clear  lakes  (Koenings  and  Edmundson  1991). 

Reporting  detection  limit  (RDL)  and  method  detection  limit  (MDL)  were 
noted  in  the  database  (Table  7-1).  The  RDL  was  the  minimum  concentration  of 
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an  analyte  that  can  be  reliably  quantified.  The  MDL  was  the  minimum 
concentration  of  an  analyte  that  can  be  detected.  Values  less  than  MDL  were 
omitted  from  the  database  and  denoted  <MDL.  Values  greater  than  MDL  but  less 
than  RDL  were  reported  to  two  significant  figures.  Values  greater  than  RDL  were 
reported  to  three  significant  figures.  Omission  of  the  <MDL  observations  created 
gaps  in  the  model-data  comparisons  and  created  the  impression  that  no 
information  was  available.  The  authors  thought  the  existence  of  observations 
below  MDL  was  important  information.  Observations  below  MDL  were  plotted 
as  MDL  in  subsequent  model-data  comparisons. 


Table  7-1 

Observed  Constituents 

Constituent 

RDL 

MDL 

Temperature,  °C 

Dissolved  Oxygen,  g  m'3 

Ammonium,  g  N  rrf3 

0.04 

0.02 

Nitrate+Nitrite,  g  N  m'3 

0.1 

0.05 

Total  Nitrogen,  g  N  rrf3 

0.1 

0.05 

Ortho  Phosphorus,  g  P  m'3 

0.005 

0.002 

Total  Phosphorus,  g  P  rrf3 

0.01 

0.005 

Chlorophyll,  mg  m'3 

0.05 

0.01 

Disk  Visibility,  m 

Fecal  Coliform,  mpn  (100  ml)'1 

Roughly  50  dissolved  and  total  organic  carbon  observations,  collected  in 
2002,  were  provided  by  the  sponsor.  Although  the  observations  were  collected 
outside  the  simulation  period,  they  were  available  for  qualitative  comparison 
with  the  model. 


Format  of  Model-Data  Comparisons 

Graphical  presentation 

Time  series.  Time  series  comparisons  of  computed  and  observed  quantities 
were  completed  for  18  stations  in  Lakes  Washington  and  Union.  Comparisons 
were  completed  at  up  to  three  depths,  depending  on  station  depth  and  sample 
protocol.  Surface,  middepth,  and  bottom  were  plotted  at  the  deepest  stations. 
Surface  and  bottom  were  plotted  at  stations  of  moderate  depth.  Surface  only  was 
plotted  for  the  shallowest  stations.  For  all  substances  except  chlorophyll,  samples 
were  compared  to  model  results  in  the  cell  corresponding  to  sample  depth. 
Chlorophyll  samples  were  compared  to  the  average  of  the  two  uppermost  model 
cells.  Depth  of  these  two  cells  approximated  the  Secchi  depth  over  which 
chlorophyll  samples  were  integrated. 

Vertical  profiles.  Vertical  profiles  were  plotted  at  three  stations  in  Lake 
Union  and  six  stations  in  Lake  Washington.  Profiles  were  prepared  on  a  monthly 
basis.  For  months  with  one  survey,  observations  were  plotted  against  model 
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mean  for  the  5 -day  period  centered  on  the  survey  date.  For  months  with  multiple 
surveys,  the  observed  mean  and  range  were  plotted  against  model  mean  for  the  5- 
day  periods  centered  on  the  survey  dates. 

Scatterplots.  Scatterplots  were  prepared  for  various  model-data 
assemblages.  In  all  cases,  observations  were  compared  to  daily-average  model 
computations  in  cells  corresponding  to  sample  depth.  Comparisons  were 
provided  for  13  individual  stations  in  Lake  Washington,  for  all  observations  in 
Lake  Washington,  for  all  observations  in  Lake  Union,  and  for  all  observations  in 
Lakes  Washington  and  Union.  An  assemblage  was  created  to  examine  potential 
differences  in  model  performance  with  depth.  Three  deep  stations  (826,  852,  890) 
were  examined  over  die  April  to  October  period  when  stratification  is  present. 
Observations  and  model  were  grouped  into  three  depth  increments: 

•  Surface  mixed  layer,  layers  1  to  6,  depth  <  1 0  m 

•  Thermocline,  layers  7  to  13,  depth  10  to  20  m 

•  Below  thermocline,  layers  14  to  41,  depth  >  20  m 

Statistical  comparisons 

Statistical  summaries  were  prepared  corresponding  to  the  scatterplot 
assemblages  for  the  following: 

a.  Mean  error. 

ME  =  ±-fi(0D-Pn)  (7-2) 

N  „=, 


in  which: 

N  =  number  of  observations 
On  =  nth  observation 

P„  -  computation  corresponding  to  nth  observation 

b.  Relative  error. 


RE  = 


|me| 

Gojljv 


(7-3) 


c.  Variance. 


VAR  = 


^  n= 1 


- ME 2 


d.  Linear  correlation. 


(7-4) 
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e.  Below  MDL.  Observations  below  the  MDL  were  treated  with  the 
following  protocol  in  computing  statistics.  If  both  observation  and 
computation  were  below  MDL,  the  observation  was  set  equal  to  the 
computation.  If  the  observation  was  below  MDL  while  the  computation 
exceeded  MDL,  the  observation  was  set  to  MDL. 


Results 

Graphical  comparisons  of  model  results  were  produced  using  Tecplot,  a 
commercially  available  graphics  package.  Results  were  plotted  one  to  a  page, 
resulting  in  several  thousand  pages  of  plots.  These  were  assembled  into 
bookmarked  PDF’s  (portable  document  format)  and  provided  to  the  project 
sponsor.  The  PDF’s  are  also  available  from  the  first  author.  Illustrations  in  this 
report  represent  a  limited,  especially  significant,  selection.  Time  series  and 
vertical  plots  are  presented  for  Station  852  (Figure  7-1).  This  centrally  located 
station  is  situated  in  the  deepest  portion  of  the  lake  and  coincides  with  the 
locations  of  studies  conducted  over  a  period  of  decades. 

Station  852 

Temperature.  Surface  temperature  at  Station  852  exhibits  an  annual  range 
of  «15  °C  (Figure  7-3).  Minimum  temperature  of  ~7  °C  occurs  in  late  February  or 
early  March.  Maximum  temperature  of  22  °C  or  higher  occurs  in  late  July  or 
early  August.  Temperature  fluctuations  are  damped,  and  the  date  of  maximum 
temperature  is  delayed  with  increasing  depth  into  the  water  column.  At  middepth 
(Figure  7-4),  the  annual  cycle  moves  through  ~5  °C  and  maximum  temperature 
occurs  in  late  October.  At  bottom  (Figure  7-5),  the  fluctuations  are  damped  to 
less  than  3  °C  and  maximum  temperature  occurs  in  late  November.  The  model 
provides  excellent  representation  of  amplitude  and  timing  of  the  annual  cycle 
although  the  temperatures  at  middepth  are  high  by  ~2  °C  in  the  latter  half  of  the 
year. 

The  lake  is  vertically  well-mixed  from  January  through  March  (Figures  7-6 
through  7-8).  Surface  warming  commences  in  April,  and  athennocline  develops 
by  June.  A  distinct  surface-mixed  layer  may  or  may  not  be  evident  in  June  and 
July.  By  August,  this  feature  is  clearly  developed.  The  surface-mixed  layer 
increases  in  depth  as  the  water  surface  cools  through  November.  By  December, 
mixing  is  nearly  complete  from  top  to  bottom.  The  model  provides  an  excellent 
representation  of  the  vertical  structure  of  the  lake.  During  the  summer  months, 
the  observations  show  a  sharper  thermocline  and  more  uniform  temperature 
below  the  thermocline  than  the  model.  The  smoother  model  profile  results  in 
higher  computed  temperatures  just  below  the  thermocline.  These  were  also 
evident  in  the  time  series  plot  at  middepth.  Extensive  experiments  were 
perfonned  with  the  hydrodynamic  model  to  tiy  to  sharpen  the  thermocline. 
Present  results  are  the  best  that  can  be  obtained.  It  is  uncertain  if  the  vertical 
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turbulence  closure  scheme  requires  revision  or  if  the  present  grid  spacing,  1.5  m, 
should  be  refined  downwards. 

Chlorophyll.  Surface  chlorophyll  at  Station  852  demonstrates  a  spring  algal 
bloom  that  commences  in  early  April  and  is  over  by  the  end  of  May  (Figure  7-9). 
During  the  bloom,  concentrations  rise  from  less  than  5  mg  m'3  to  more  than  15 
mg  m'3.  Following  the  bloom,  concentrations  return  to  less  than  5  mg  m'  . 
Observations  suggest  the  occurrence  of  an  irregular,  lesser  bloom  in  late  October 
to  early  November.  The  model  performs  well  in  representing  the  timing  and 
magnitude  of  the  spring  bloom  as  well  as  the  occurrence  and  magnitude  of  the 
fall  bloom. 

The  predation  term  was  found  to  be  the  key  to  representing  the  termination  of 
the  spring  bloom.  A  piece-wise  first-order  predation  term  was  used  in  which 
predation  increased  from  0.045  day"1  (November-April)  to  0.22  day  1  (May- 
October).  Employment  of  the  lower  rate  throughout  the  year  resulted  in  a  spring 
bloom  of  longer  duration  and  in  summer  chlorophyll  concentrations  roughly 
double  observed  concentrations  (Figure  7-10).  Computations  suggest  the 
influence  of  top-down  control  of  phytoplankton  during  late  spring  and  summer. 
Explicit  representation  of  zooplankton  will  be  a  worthwhile  improvement  in  the 
model. 

Disk  visibility.  The  model  computes  the  light  attenuation  coefficient  via  an 
additive  relationship: 

Ke  =  a  +  b  ■  ISS  +  c  ■  Chi  (7-6) 


in  which: 

a  =  background  light  attenuation  (0.29  m"1) 

b  =  coefficient  that  relates  extinction  to  fixed  solids  concentration 
(0.08  m2  g"1) 

ISS  =  fixed  (inorganic)  solids  (g  m"3) 

c  =  coefficient  that  relates  extinction  to  chlorophyll  concentration 
(0.02  m2  mg"1) 

Chi  =  chlorophyll  concentration  (mg  m"3) 

Coefficients  a  and  b  were  assigned  values  based  on  experience  in  other 
systems  (Cerco  and  Noel  in  preparation).  Coefficient  c  was  assigned  based  on  a 
summary  of  observed  values  (Pennock  1985).  Attenuation  by  organic  solids  is 
incorporated  in  the  chlorophyll  attenuation  coefficient  since  organic  solids  are 
largely  associated  with  phytoplankton.  Computed  extinction  is  converted  to  disk 
visibility,  for  comparison  with  observations,  via  Equation  7-1. 

Background  attenuation,  attributed  primarily  to  color,  is  a  major,  consistent 
component  of  computed  attenuation  (Figure  7-11).  In  the  absence  of  any  solids, 
disk  visibility  in  Lake  Washington  would  be  ~6  m.  Computed  fixed  solids  vary 
from  0  to  ~1  g  m"3,  with  highest  concentrations  occurring  in  late  winter  and  early 
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spring,  coincident  with  maximum  runoff  (Figure  7-12).  Fixed  solids  contribute 
=0.08  m1  to  computed  light  attenuation  during  their  period  of  maximum 
concentration.  They  reduce  disk  visibility  by  1 .5  m,  at  most.  During  the  spring 
bloom,  attenuation  from  chlorophyll  rivals  color  as  the  largest  contributor  to  the 
total  (Figure  7-1 1).  Chlorophyll  is  the  major  contributor  to  variance  in  light 
attenuation  and,  hence,  disk  visibility.  The  range  in  computed  attenuation  is  from 
nearly  zero  to  0.3  m1.  During  the  bloom,  chlorophyll  reduces  computed  disk 
visibility  by  =2.5  m. 

Computed  disk  visibility  (Figure  7-13)  is  virtually  the  inverse  of  computed 
chlorophyll.  Maximum  computed  visibility  occurs  in  late  autumn,  when 
chlorophyll  is  low  and  before  winter  runoff  events.  Runoff-induced  solids 
diminish  visibility  by  =1  m  from  November  through  February.  Disk  visibility 
reaches  its  minimum,  =3  m,  in  late  March  to  early  April,  during  the  spring 
bloom,  and  hovers  around  5  m  during  the  summer  months.  The  computed  pattern 
agrees  reasonably  well  with  the  observations,  which  show  a  consistent  minimum 
coincident  with  die  spring  bloom  and  suggest  a  maximum  in  late  autumn.  The 
observations  show  a  great  deal  of  scatter  about  the  model  value  and  a  great  deal 
of  scatter  in  general.  The  loose  agreement  between  model  and  observed  visibility 
is  attributed  to  three  factors:  the  inherent  uncertainty  in  the  disk  visibility 
measure,  lack  of  solids  data  for  calibration  of  this  attenuating  factor,  and  the 
conversion  from  attenuation  to  disk  visibility.  Aside  from  uncertainty  in  the 
Poole  and  Atkins  relationship,  the  inverse  relationship  transforms  small 
variations  in  computed  attenuation  to  large  variations  in  disk  visibility.  The 
model  is  believed  to  properly  represent  the  magnitude  of  disk  visibility  and  the 
processes  that  contribute  to  seasonal  variations  in  this  observed  quantity. 

Nitrogen.  The  preponderance  of  ammonium  observations  at  Station  852  is 
below  MDL  of  0.02  g  m’3  (Figures  7-14  through  7-16).  Highest  concentrations 
are  roughly  double  the  MDL.  The  available  observations  reveal  little  regular 
temporal  or  vertical  variability.  Certainly,  the  model  reflects  the  magnitude  and 
range  of  the  observations.  The  observations  for  the  latter  2  years  of  the 
simulation  hint  at  higher  concentrations  in  spring,  perhaps  indicating 
zooplankton  excretion  following  predation  on  the  spring  bloom.  The  model 
suggests  a  later  peak,  originating  in  predation  on  the  summer  algal  population. 

Although  the  nitrate  observations  at  this  station  are  sparse,  the  available 
observations  here  and  elsewhere  indicate  clear  temporal  and  spatial  trends.  At  the 
surface  of  Station  852  (Figure  7-17),  peak  nitrate  concentrations,  =0.3  g  m'3, 
occur  in  late  winter  and  early  spring.  In  these  months,  nitrate  supplied  by  runoff 
lacks  major  internal  sinks.  As  the  water  warms  and  days  get  longer,  surface 
nitrate  is  taken  up  by  phytoplankton  until  a  minimum,  0.02  g  m'3  or  less,  is 
reached  in  late  summer.  Diminished  uptake  in  autumn  coupled  with  increasing 
runoff  returns  surface  nitrate  to  its  maximum  concentration  in  late  winter.  The 
model  performs  admirably  in  reproducing  the  magnitude  and  timing  of  the 
annual  cycle. 

Computed  nitrate  cycling  at  middepth  is  influenced  by  algal  uptake  and  by 
the  vertical  physical  structure  of  the  lake  (Figure  7-18).  A  minimum  occurs  in 
spring,  coincident  with  the  spring  bloom  (Figure  7-19).  Thereafter,  nitrate 
concentration  increases,  driven  by  sequential  hydrolysis  of  organic  nitrogen  and 
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nitrification  of  ammonium  and,  perhaps,  by  diffusion  from  below.  A  second 
minimum  occurs  in  late  autumn,  following  fall  overturn,  when  nitrate-depleted 
surface  water  is  mixed  downwards.  Nitrate  again  increases,  due  to  runoff,  until 
the  recurrence  of  the  spring  bloom.  Available  data  indicate  the  model  is 
performing  well  in  reproducing  the  annual  cycle  at  middepth,  especially  the  rapid 
decrease  in  nitrate  following  fall  overturn. 

Computed  nitrate  concentration  near  the  bottom  (Figure  7-20)  hints  at  an 
annual  cycle,  but  the  phenomenon  is  difficult  to  verify.  The  model  profile  (Figure 
7-18)  indicates  that  nitrate  is  lower  immediately  above  the  bottom  than  near 
middepth  from  May  through  December.  This  shape  is  sometimes  apparent  in  the 
observations  but  not  as  pronounced  as  in  the  model.  Overall,  computed 
concentrations  in  the  model  are  lower  than  the  available  observations. 
Observations  range  from  ~0.3  to  0.4  g  m'3,  while  comparable  computations  are  as 
much  as  0.2  g  m'3  lower. 

Total  nitrogen  observed  at  the  surface  of  Station  852  follows  an  annual 
pattern  that  closely  corresponds  to  nitrate,  which  comprises  a  major  fraction  of 
the  total  (Figure  7-21).  Computed  cycling  matches  the  observations  well 
although  the  model  oscillations  are  damped  relative  to  the  observations. 

At  middepth,  observed  and  computed  total  nitrogen  demonstrate  similar 
behavior  and  magnitude  (Figure  7-22).  The  total  dips  to  an  annual  minimum,  0.2 
to  0.3  g  m'3,  following  the  mixing  of  nitrate-depleted  water  downward  during  fall 
turnover.  By  late  winter,  total  nitrogen  recovers  to  a  value  that  remains  nearly 
constant  until  the  next  turnover. 

At  bottom,  observed  total  nitrogen  varies  from  0.3  to  0.6  g  m'3  with  no 
obvious  repeated  cycle  (Figure  7-23).  The  model  demonstrates  weak  cycling  with 
a  peak  in  midsummer.  The  peak  is  composed  of  ~0. 1  g  m‘3  particulate  nitrogen 
that  settles  to  the  bottom  following  the  spring  bloom. 

Phosphorus.  At  the  surface  of  Station  852,  ortho-phosphorus  shows  annual 
cycling  (Figure  7-24)  very  like  nitrate.  Highest  concentrations  occur  in  late 
winter  to  early  spring  when  runoff  is  high  and  internal  demand  is  low.  Ortho¬ 
phosphorus  is  almost  instantly  depleted  by  the  spring  algal  bloom  and  remains 
low  throughout  the  summer.  In  autumn,  in  the  absence  of  algal  uptake,  ortho¬ 
phosphorus  recovers  to  its  annual  maximum.  Observed  concentrations  range 
from  below  MDL  (0.002  g  m'3)  to  -0.02  g  m‘3.  The  model  performs  well  in 
reproducing  the  observed  magnitude  and  cycling. 

At  middepth,  regular  cycling  of  ortho-phosphorus  is  difficult  to  detect  in  the 
observations  although  lower  values  usually  occur  in  late  summer  versus  higher 
values  in  late  winter  to  early  spring  (Figure  7-25).  The  observed  range  varies 
from  below  MDL  (0.002  g  m'3)  to  -0.02  g  m"3.  The  model  exhibits  cycling  much 
like  nitrate  at  the  same  depth.  One  minimum  occurs  following  uptake  by  the 
spring  bloom.  Thereafter,  ortho-phosphorus  recovers  to  an  annual  maximum.  A 
second  minimum  occurs  when  phosphorus-depleted  surface  water  is  mixed 
downward  during  fall  overturn.  Runoff  and  absence  of  uptake  create  a  secondary 
maximum  prior  to  the  demand  created  by  the  spring  bloom. 
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Although  the  model  cycle  is  readily  explainable,  it  is  not  apparent  in  the  data. 
In  particular,  the  recover,'  following  the  spring  bloom  is  missing  from  the 
observations.  This  phenomenon  of  excess  ortho-phosphorus  is  evident  in  other 
model  applications  (e.g.,  Cerco  and  Noel  in  preparation).  The  problem  originates 
in  two  features  of  conventional  phosphorus  models. 

In  the  conventional  modeling  approach,  dissolved  organic  phosphorus 
mineralizes  to  ortho-phosphorus  at  a  first-order  rate  regardless  of  demand.  No 
mechanistic  explanation  is  provided  to  explain  why  a  fixed  fraction  of  organic 
phosphorus  is  mineralized  on  a  daily  basis.  Conventional  models  include  no  sink 
for  ortho-phosphorus  except  algal  uptake.  Ortho-phosphorus  produced  in  the 
dark  by  first-order  decay  lasts  eternally. 

This  model  addresses  the  issue  of  excess  ortho-phosphorus  production  by 
relating  mineralization  rate  to  algal  biomass  and  nutrient  deficit.  Mineralization  is 
reduced  when  biomass  is  low7  and/or  when  ortho-phosphorus  is  abundant  relative 
to  algal  requirements  (see  Chapter  3).  In  an  earlier  approach,  DiToro  and  Matysik 
(1979)  also  related  mineralization  to  algal  biomass.  While  the  supply  side  can  be 
addressed  by  empirical  formulations,  demand  by  biota  other  than  phytoplankton 
seems  to  be  universally  absent  from  eutrophication  models.  A  cursoiy  literature 
search  indicates  that  heterotrophic  bacteria  compete  with  phytoplankton  for 
phosphate  (Thingstad  et  al.  1993)  and  may  draw  phosphate  down  to  limiting 
levels  (Zweifel  et  al.  1993).  While  an  explicit  representation  of  bacterial 
phosphorus  uptake  is  a  major  undertaking,  an  empirical  approach  should  be 
feasible.  One  approach  is  to  relate  phosphorus  uptake  to  organic  carbon 
respiration,  which  is  a  bacterial  process. 

On  the  bottom,  computed  ortho-phosphorus  demonstrates  an  annual  cycle 
that  is  representative  of  the  observations  (Figure  7-26).  The  annual  minimum 
occurs  following  the  autumn  turnover  when  phosphorus-depleted  water  is  mixed 
downward.  Thereafter,  ortho-phosphorus  recovers  to  an  annual  maximum  prior 
to  the  next  turnover.  The  cycling  and  range,  0.01  to  0.03  g  m'3,  of  the 
observations  are  well  represented  by  the  model. 

Vertical  profiles  for  ortho-phosphorus  (Figures  7-27  through  7-29)  provide  a 
different  perspective  but  indicate  the  same  behaviors  as  the  time  series.  At  the 
surface,  the  model  represents  the  depletion  of  ortho-phosphorus  during  the  late 
spring  and  summer  months.  Below  the  tbermocline,  however,  ortho-phosphorus 
is  trapped  and  forms  a  bulge  that  persists  until  turnover  in  January  or  February. 

Observed  and  computed  total  phosphorus  levels  in  surface  waters  follow  a 
seasonal  cycle  (Figure  7-30).  Maximum  concentrations,  «0.03  g  m"3,  occur  in  late 
winter,  when  runoff-induced  loads  are  largest.  Diminished  loads  and  loss  through 
sinking  of  algal  detritus  produce  a  minimum,  ~0.005  g  m'3,  in  late  summer. 

Seasonal  cycling  is  difficult  to  perceive  in  the  observations  at  the  middepth 
(Figure  7-3 1)  and  near  bottom  (Figure  7-32).  More  distinctive  cycling  is  apparent 
in  the  model.  Near-bottom,  maximum  computed  concentrations,  -0.05  g  m'3, 
occur  in  late  spring,  following  deposition  of  the  algal  bloom.  Thereafter 
concentration  tails  off  slowly  until  the  winter  overturn.  At  die  overturn, 
computed  total  phosphorus  declines  sharply  as  phosphorus-depleted  surface 
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water  is  mixed  with  bottom  water.  The  annual  minimum,  -0.02  g  m'3,  occurs  in 
late  winter  and  persists  until  deposition  of  the  spring  algal  bloom.  Cycling  of 
computed  total  phosphorus  at  middepth  is  damped  relative  to  bottom  water.  The 
predominant  feature  is  a  late-fall  dip  of —0.01  g  m 3  as  surface  water  is  mixed 
downwards. 

In  the  vertical  dimension  (Figures  7-33  through  7-35),  cycling  between  well- 
mixed  and  stratified  conditions  is  evident  in  the  observations.  Total  phosphorus 
is  well  mixed  from  top  to  bottom  from  January  through  March  or  April. 
Thereafter,  concentrations  near  the  bottom  increase  to  the  highest  values 
observed  at  the  station,  occasionally  in  excess  of  0.05  g  m 3  from  August  through 
November.  Comparison  with  the  ortho-phosphorus  profiles  indicates  that  the 
preponderance  of  the  total  is  ortho-phosphorus,  most  likely  released  from  bottom 
sediments. 

Organic  carbon.  Organic  carbon  observations  were  not  collected  during  the 
model  application  period.  A  limited  number  of  observations  were  available  from 
2002.  For  qualitative  comparison  with  the  model,  these  observations  were  plotted 
on  dates  in  1995-1997  corresponding  to  the  month  and  day  of  collection  in  2002. 

The  observations  suggest  particulate  organic  carbon  near  the  surface  peaks  at 
more  than  1  g  m'3  during  the  spring  bloom  (Figure  7-36).  Otherwise 
concentrations  are  —0.5  g  m 3.  Model  computations  reproduce  this  behavior.  Near 
the  bottom,  observed  particulate  organic  carbon  is  «0.5  g  m'3  or  less  (Figure  7- 
37).  The  observations  hint  at  a  late  winter  peak,  but  absence  of  a  similar  peak  in 
observations  from  other  stations  suggests  that  the  apparent  pattern  at  Station  852 
cannot  be  discerned  from  random  variance.  The  model  shows  a  definite  seasonal 
signal.  Near-bottom  particulate  organic  carbon  peaks  in  late  spring  when  detritus 
from  the  algal  bloom  settles  to  the  bottom.  The  computed  maximum,  more  than 
1.5  g  m'3,  exceeds  values  in  the  limited  observations.  Otherwise  the  computations 
reflect  the  observed  magnitude  and  range. 

Observations  indicate  surface  dissolved  organic  carbon  ranges  from  3  to  4  g 
m'3  with  the  highest  concentration  occurring  in  late  summer  (Figure  7-38).  The 
late  summer  maximum  at  Station  852  is  consistent  with  a  maximum  at  other 
stations.  The  model  also  shows  a  summer  maximum,  associated  with  algal 
excretion  and  with  byproducts  of  predation.  The  computed  maximum  is  in  the 
observed  range.  The  computed  minimum,  — 2  g  m 3,  is  less  than  observed,  but  the 
sparse  observations  preclude  any  conclusions  about  model  behavior  in  this 
regard. 

Near-bottom  observed  dissolved  organic  carbon  is  ~3  g  m'3  (Figure  7-39). 
Computed  dissolved  organic  carbon  is  1  to  2  g  m 3  less  than  observed.  The  lower 
computed  values  are  consistent  with  higher  computed  respiration,  as  will  be  seen 
in  the  next  chapter.  Cyclic  behavior  is  apparent  in  the  model  near  bottom. 

Highest  concentrations  occur  in  winter  when  dissolved  organic  carbon  in  runoff 
is  mixed  downward.  Thereafter,  dissolved  organic  carbon  slowly  decays  until  the 
subsequent  turnover  event. 

Dissolved  oxygen.  Observed  and  computed  dissolved  oxygen  levels  at 
Station  852  show  strong  seasonal  cycling  at  all  depths.  Near  the  surface, 
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observed  dissolved  oxygen  reaches  13  to  14  g  m'3  in  late  spring  (Figure  7-40). 
These  concentrations  are  supersaturated,  indicating  production  by  the  spring 
bloom  (Figure  7-41).  Dissolved  oxygen  declines  but  remains  supersaturated 
through  late  summer  when  the  annual  minimum  of  8  to  9  g  m'3  is  attained. 
Thereafter,  dissolved  oxygen  increases,  but  remains  subsaturated  until  the 
subsequent  spring  algal  bloom.  The  model  performs  well  in  reproducing  the 
observed  concentrations  as  well  as  the  cycling  from  supersaturated  to 
subsaturated  conditions  (Figure  7-42). 

At  middepth,  observed  and  computed  dissolved  oxygen  levels  are  1  to 
2  g  m'3  less  than  in  surface  waters  (Figure  7-43).  The  decline  from  the  spring 
maximum  is  smoother  than  in  surface  waters,  indicative  of  respiratory  activity 
rather  than  a  sharp  decline  in  algal  oxygen  production.  The  annual  minimum, 

~7  g  m'3,  is  short-lived;  dissolved  oxygen  increases  rapidly  when  the  surface- 
mixed  layer  moves  downward  in  autumn. 

The  largest  annual  range  in  dissolved  oxygen  occurs  near  the  bottom  (Figure 
7-44).  Computed  and  observed  annual  maximums  greater  than  1 1  g  m'3  occur  in 
spring,  coincident  with  the  algal  bloom  and  prior  to  stratification.  As  the  lake 
stratifies,  dissolved  oxygen  declines  to  ~4  g  m'3  in  late  autumn,  then  recovers 
rapidly  when  the  lake  completely  overturns  in  winter. 

The  most  interesting  feature  of  the  vertical  dissolved  oxygen  structure  is  the 
metalimnetic  oxygen  minimum  that  occurs  from  July  through  October  or 
November  (Figures  7-45  through  7-47).  The  processes  that  create  the  minimum 
are  varied  (Wetzel  1975).  Two  processes  that  have  been  assigned  roles  in  the 
fonnation  of  the  minimum  are  absent  in  this  model:  density-induced  differential 
settling  of  particulate  organic  matter  (Wetzel  1975)  and  copepod  respiration 
(Shapiro  1960).  The  modeled  minimum  most  likely  occurs  at  the  location  where 
algal  oxygen  production  is  absent  yet  temperature  is  sufficient  to  produce  high 
respiration  of  dissolved  organic  carbon.  At  lower  depths,  respiration  is  reduced 
by  lower  temperatures. 

Fecal  coliform.  Observed  fecal  coliform  levels  in  surface  waters  range  up  to 
60  mpn  (100  ml)'1  (Figure  7-48).  The  preponderance  of  values,  however,  is  less 
than  20  mpn  (100  ml)'1.  Minimum  observed  values  occur  in  the  summer  and 
early  autumn,  most  likely  because  storm  events  are  absent  and  temperature- 
induced  mortality  is  high.  In  late  autumn  and  winter,  fecal  coliform 
concentrations  are  randomly  distributed.  The  maximums  are,  no  doubt, 
associated  with  loading  events.  The  model  demonstrates  a  seasonal  behavior 
consistent  with  the  observations  although  one-to-one  correspondence  with 
individual  observations  is  absent.  Lowest  modeled  concentrations  occur  in 
summer.  As  surface  waters  cool,  mortality  is  reduced  and  concentrations  increase 
until  winter  overturn.  The  overturn  causes  a  dramatic  decrease  as  surface  waters 
are  mixed  downwards.  Thereafter,  computed  concentrations  increase  erratically 
until  warm  summer  temperatures  and  absence  of  runoff  events  create  the  summer 
minimum  concentrations. 

Modeled  concentrations  at  this  location  are  less  than  20  mpn  (100  ml)'1,  in 
agreement  with  the  preponderance  of  the  observations.  The  highest  observed 
values  are  not  represented  in  the  model.  The  computations  are  affected  by  the 
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crude  loading  information;  e.g.,  all  combined-sewer  overflow  events  are 
characterized  by  a  single  concentration.  Improved  calculations  of  fecal  coliform 
will  require  improved  loading  estimates. 

Statistics  and  scatterplots 

Statistical  summaries  of  model  results  can  provide  insights  not  evident  in 
graphical  comparisons  of  observed  and  computed  properties.  Model  statistics 
present  problems,  as  well.  No  standard  set  of  statistics  exists.  Neither  are  there 
standard  criteria  for  acceptance  of  model  results.  In  the  absence  of  formal 
hypothesis  testing,  judgment  must  be  employed  in  determining  whether 
differences  in  statistics  for  various  model  runs,  substances,  or  groupings  are 
“significant”  or  not. 

Comprehensive  model  statistics  are  presented  in  Table  7-2.  Two  very  useful 
statistics  are  mean  error  (ME)  and  relative  error  (RE).  Mean  error  indicates 
whether  model  results  are  higher  or  lower  than  observations,  on  average.  Relative 
error  is  useful  for  comparisons  between  various  substances  or  model  applications 
to  various  systems.  This  discussion  focuses  on  these  statistics  as  applied  to  Lake 
Washington. 


Table  7-2 

Summary  Statistics 

Station 

Constituent 

Relative  Error 

Mean  Error 

Correlation 

Variance 

Station  852 

Ammonium  (mg/L) 

0.002 

-0.00005 

0.33 

0.00008 

Lake  Washington  Stations 

0.026 

0.0006 

0.21 

0.00013 

Lake  Union  Stations 

0.266 

0.00855 

0.17 

0.00052 

826,  859,  890  (Apr-Oct)  Level  1-6 

0.026 

-0.00061 

0.25 

0.00013 

826,  859,  890  (Apr-Oct)  Level  7-13 

0.019 

0.00043 

0.29 

0.00012 

826,  859,  890  (Apr-Oct)  Level  14-41 

0.096 

0.00195 

0.18 

0.00011 

Station  852 

Chlorophyll  (ugm/L) 

0.341 

0.72 

0.80 

5.5 

Lake  Washington  Stations 

0.047 

0.20 

0.70 

11.9 

Lake  Union  Stations 

0.123 

0.59 

0.58 

15.5 

826,  859,  890  (Apr-Oct)  Level  1-6 

0.136 

0.65 

0.80 

10.4 

Station  852 

Dissolved  Oxygen  (mg/L) 

0.038 

-0.34 

0.88 

0.91 

Lake  Washington  Stations 

0.049 

-0.46 

0.85 

1.06 

Lake  Union  Stations 

0.124 

-1.14 

0.82 

2.58 

Surface  <  10m 

0.059 

-0.54 

0.88 

0.79 

Thermocline 

0.074 

-0.66 

0.88 

0.76 

Deep  Wkter  >  20m  ! 

0.041 

-0.36 

0.88 

0.79 

Station  852 

Nitrate- Nitrite  (mg/L) 

0.118 

0.03006 

0.87 

0.00255 

Lake  Washington  Stations 

0.120 

0.02238 

0.88 

0.00266 

Lake  Union  Stations 

0.212 

0.02738 

0.92 

0.00153 

826,  859,  890  (Apr-Oct)  Level  1-6 

0.131 

0.02658 

0.90 

0.00206 

826,  859,  890  (Apr-Oct)  Level  7-13 

0.134 

0.03164 

0.86 

0.00201 

826,  859,  890  (Apr-Oct)  Level  14-41 

0.151 

0.04287 

0.65 

0.00171 

(Continued)  | 
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Table  7-2  (Concluded) 

Station 

Constituent 

Relative  Error 

Mean  Error 

Correlation 

Variance 

Station  852 

Ortho-Phosphate  (mg/L) 

0.306 

-0.00308 

0.67 

0.00004 

Lake  Washington  Stations 

0.197 

-0.00167 

0.69 

0.00003 

Lake  Union  Stations 

0.326 

0.00246 

0.33 

0.00008 

826,  859,  890  (Apr-Oct)  Level  1-6 

0.374 

-0.00311 

0.66 

0.00004 

826,  859,  890  (Apr-Oct)  Level  7-1 3 

0.466 

-0.00452 

0.62 

0.00004 

826,  859,  890  (Apr-Oct)  Level  14-41 

0.331 

-0.00455 

0.43 

0.00005 

Station  852 

Temperature  (Celcius) 

0.071 

-0.80 

0.97 

1.23 

Lake  Washington  Stations 

0.089 

-1.09 

0.96 

1.70 

Lake  Union  Stations 

0.034 

-0.49 

0.99 

0.97 

Surface  <  10m 

0.094 

-1.13 

0.96 

1.63 

Thermocline 

0.111 

-1 .22 

0.93 

1.73 

Deep  Wkter  >  20m 

0.116 

-0.94 

0.93 

0.64 

Station  852 

Total  Nitrogen  (mg/L) 

0.097 

0.03513 

0.69 

0.00339 

Lake  Washington  Stations 

0.071 

0.02528 

0.63 

0.00705 

Lake  Union  Stations 

0.153 

0.05272 

0.68 

0.00435 

826,  859,  890  (Apr-Oct)  Level  1-6 

0.061 

0.02168 

0.72 

0.00391 

826,  859,  890  (Apr-Oct)  Level  7-13 

0.077 

0.02912 

0.64 

0.00318 

826,  859,  890  (Apr-Oct)  Level  14-41 

0.101 

0.04085 

0.26 

0.00303 

Station  852 

Total  Phosphorus  (mg/L) 

0.044 

0.00087 

0.49 

0.00012 

Lake  Washington  Stations 

0.164 

0.00345 

0.35 

0.00018 

Lake  Union  Stations 

0.398 

0.00842 

0.28 

0.0001 

826,  859,  890  (Apr-Oct)  Level  1-6 

0.051 

0.00097 

0.40 

0.00013 

826,  859,  890  (Apr-Oct)  Level  7-13 

0.050 

-0.001 

0.45 

0.00012 

826,  859,  890  (Apr-Oct)  Level  14-41 

0.069 

-0.001 61 

0.41 

0.00014 

Temperature.  Overall,  computed  temperature  exceeds  observed  by  roughly 
1  °C.  The  excess  is  a  little  larger  at  the  surface  and  in  the  thermocline  than  in 
deep  water.  The  scatterplot  (Figure  7-49)  indicates  the  largest  excess 
temperatures  occur  when  observed  temperatures  exceed  10  °C.  Relative  error  in 
temperature  computations  is  less  than  9  percent. 

Chlorophyll.  Overall,  observed  chlorophyll  exceeds  computed  by  less  than 
0.2  mg  m'3.  Insufficient  data  exist  to  examine  performance  as  a  function  of  depth. 
The  scatterplot  (Figure  7-50)  indicates  that  magnitude  of  computed  error 
increases  as  a  function  of  magnitude  of  observed  chlorophyll.  Similar  behavior 
has  been  demonstrated  in  multiple  alternate  model  applications.  Relative  error  in 
chlorophyll  computations  is  ~5  percent. 

Nitrogen.  Observed  total  nitrogen  exceeds  computed  by  0.025  g  m'3,  on 
average.  The  computed  shortfall  in  deep  water  is  nearly  double  the  shortfall  at  the 
surface.  Computed  nitrate  is  0.022  g  m"3  less  than  observed,  on  average,  with 
larger  discrepancies  in  the  deep  waters.  Effectively,  the  entire  deficit  in  total 
nitrogen  is  in  the  nitrate  fraction.  Observed  ammonium  exceeds  computed  by 
0.0006  g  m'3,  on  average.  Consistent  with  the  other  components,  mean  error  is 
greater  in  deep  waters.  Relative  errors  for  total  nitrogen,  nitrate,  and  ammonium 
are  7  percent,  12  percent,  and  3  percent,  respectively.  Scatterplots  (Figures  7-51 
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through  7-53)  indicate  that  errors  in  total  nitrogen  and  its  components  are 
homoscedastic.  That  is,  the  magnitude  of  the  error  does  not  increase  in  proportion 
to  the  magnitude  of  the  observation. 

Phosphorus,  Observed  total  phosphorus  exceeds  computed  by  0.0035  g  m'3, 
on  average.  The  analyses  of  mean  error  by  depth  indicates  an  opposing  trend. 
Computed  total  phosphorus  exceeds  observed  by  up  to  0.0016  g  m  in  deep 
waters.  The  analyses  by  depth  were  restricted  to  the  April -October  period.  The 
opposing  mean  errors  indicate  that  total  phosphorus  is  undercomputed  in  winter 
but  overcomputed  in  deep  water  during  summer.  Relative  error,  for  all 
observations,  is  16  percent.  The  scatterplot  (Figure  7-54)  indicates  errors  are 
largely  homoscedastic,  although  the  magnitude  of  some  extreme  observations  is 
not  attained  by  the  model  at  all. 

Computed  ortho-phosphorus  exceeds  observed,  on  average  by  0.0017  g  m'3. 

A  small  tendency  exists  for  errors  to  be  larger  in  the  thermocline  and  deep  water, 
and  in  summer  rather  than  winter.  These  properties  reflect  the  trapping  of  ortho¬ 
phosphorus  at  middepth  during  summer  stratification.  Relative  error  is  nearly 
20  percent.  The  scatterplot  (Figure  7-55)  indicates  no  apparent  relationship 
between  error  and  magnitude  of  observations. 

Dissolved  oxygen.  Computed  dissolved  oxygen  exceeds  observed  by  less 
than  0.5  g  m'3,  on  average.  Largest  discrepancies  appear  to  be  in  the  thermocline 
region  although  it  is  difficult  to  discern  if  this  apparent  property  is  significant.  No 
apparent  relationship  exists  between  error  and  concentration  (Figure  7-56). 
Relative  error  is  5  percent. 

Comparison  with  other  systems 

Relative  error  was  compared  for  a  few  key  parameters  with  relative  error  in 
two  recent  model  applications,  Chesapeake  Bay  (Cerco  and  Noel  in  preparation) 
and  the  lower  St.  Johns  River,  Florida  (Tillman  et  al.  in  preparation).  Chesapeake 
Bay  is  a  partially  mixed  estuary  located  in  the  mid-Atlantic  region.  The  lower 
St.  Johns  is  a  partially  to  well-mixed  estuary  located  in  northeast  Florida. 

Relative  error  for  Lake  Washington  (Figure  7-57)  is  much  less  than  for  the 
estuarine  applications.  The  estuaries  demonstrate  highest  relative  error  for 
chlorophyll  while  Lake  Washington  demonstrates  highest  relative  error  for 
phosphorus. 

The  lower  error  for  Lake  Washington  is  suspected  to  reflect  the  properties  of 
a  lake  relative  to  an  estuary.  The  lake  is  a  basin.  Material  introduced  to  the  lake 
must  either  settle  or  flow  out  over  the  dam.  An  estuary  is  open  at  the  lower  end. 
Tides  move  material  upstream  and  downstream  twice  a  day.  Net  transport  across 
the  open  boundary  is  especially  difficult  to  quantify. 

Computation  of  phosphorus  in  the  estuaries,  especially  Chesapeake  Bay,  has 
been  found  to  be  equally  vexing  as  in  Lake  Washington.  The  nature  of  the 
problem  is  the  same.  Without  bacterial  uptake,  dissolved  phosphate  accumulates 
to  excess  levels.  In  the  estuaries,  however,  chlorophyll  relative  error  exceeds 
phosphorus  relative  error.  The  reason  for  the  lower  chlorophyll  error  in  Lake 
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Washington  is  not  apparent.  Both  estuaries  experience  spring  algal  blooms  as 
does  Lake  Washington.  Consequently,  similar  relative  errors  would  be  expected. 

The  dissolved  oxygen  relative  error  requires  some  explanation.  Chesapeake 
Bay  experiences  bottom-water  anoxia  during  summer  months.  Since  prediction 
of  this  anoxia  is  critical,  the  statistical  computation  is  restricted  to  summer, 
bottom  dissolved  oxygen.  Statistics  for  the  St.  Johns  and  Lake  Washington 
include  all  observations.  These  tend  to  be  improved,  relative  to  Chesapeake  Bay, 
since  they  include  surface  and  cold-weather  samples  in  which  oxygen  tends  to 
equilibrate  with  the  atmosphere.  Computation  of  dissolved  oxygen  under  these 
circumstances  is  largely  a  function  of  the  algorithm  used  to  compute  saturation 
dissolved  oxygen  concentration.  Also,  the  denominator  in  the  relative  error 
computation,  essentially  saturation  concentration,  is  larger  than  the  mean  bottom 
dissolved  oxygen  that  forms  the  denominator  in  the  Chesapeake  Bay  statistics. 
Still,  relative  dissolved  oxygen  error  is  lower  in  Lake  Washington  than  in  the 
St.  Johns,  for  which  the  statistic  includes  the  entire  population  of  dissolved 
oxygen  observations. 

Virtually  all  CE-QUAL-ICM  applications  to  date  have  been  to  estuaries  and 
coastal  lagoons.  Lake  Washington  is  the  first  lake.  While  it  cannot  be  generalized 
that  lake  applications  have  lower  errors  than  estuarine  applications,  this  study 
compares  favorably  with  the  best  alternate  applications. 
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Figure  7-1 .  Lake  Washington  sample  stations 
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Figure  7-2.  Lake  Union  sample  stations 
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Figure  7-3.  Observed  and  computed  surface  temperatures  at  Station  852 
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Figure  7-4.  Observed  and  computed  middepth  temperatures  at  Station  852 
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Figure  7-5.  Observed  and  computed  bottom  temperatures  at  Station  852 
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Figure  7-6.  Observed  and  computed  temperature  profiles  at  Station  852  for  1995 
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Figure  7-7.  Observed  and  computed  temperature  profiles  at  Station  852  for  1996 
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Figure  7-8.  Observed  and  computed  temperature  profiles  at  Station  852  for  1997 
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Figure  7-13.  Observed  and  computed  Secchi  depth  at  Station  852 


Figure  7-14.  Observed  and  computed  surface  ammonium  at  Station  852 
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Figure  7-15.  Observed  and  computed  middepth  ammonium  at  Station  852 


Figure  7-16.  Observed  and  computed  bottom  ammonium  at  Station  852 
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Figure  7-18.  Observed  and  computed  nitrate  profiles  at  Station  852  for  1997 
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ire  7-21 .  Observed  and  computed  surface  total  nitrogen  at  Station  852 
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7-23.  Observed  and  computed  bottom  total  nitrogen  at  Station  852 
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Figure  7-24.  Observed  and  computed  surface  ortho-phosphorus  at  Station  852 
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Figure  7-25.  Observed  and  computed  middepth  ortho-phosphorus  at  Station  852 
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Figure  7-28.  Observed  and  computed  ortho-phosphorus  profiles  at  Station  852  for  1996 
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Figure  7-29.  Observed  and  computed  ortho-phosphorus  profiles  at  Station  852  for  1997 
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Figure  7-30.  Observed  and  computed  surface  total  phosphorus  at  Station  852 
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Figure  7-31.  Observed  and  computed  middepth  total  phosphorus  at  Station  852 
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Figure  7-33.  Observed  and  computed  total  phosphorus  profiles  at  Station  852  for  1995 
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Figure  7-34.  Observed  and  computed  total  phosphorus  profiles  at  Station  852  for  1996 
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Figure  7-35.  Observed  and  computed  total  phosphorus  profiles  at  Station  852  for  1997 
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ire  7-38.  Observed  and  computed  surface  dissolved  organic  carbon  at  Station 
852 


LAKE  WASHINGTON  (jul23  2nd) 
Dissolved  Organic  Carbon 
Station  852  Bottom 


Figure  7-39.  Observed  and  computed  bottom  dissolved  organic  carbon  at  Station 
852 
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Figure  7-40.  Observed  and  computed  surface  dissolved  oxygen  at  Station  852 


Figure  7-41 .  Observed  and  saturated  surface  dissolved  oxygen  concentrations  at 
Station  852 
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re  7-42.  Modeled  and  saturated  surface  dissolved  oxygen  concentrations  at 
Station  852 
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Figure  7-43.  Observed  and  computed  middepth  dissolved  oxygen  at  Station  852 
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LAKE  WASHINGTON  (jul23  2nd) 
Dissolved  Oxygen 
Station  852  Bottom 


Figure  7-44.  Observed  and  computed  bottom  dissolved  oxygen  at  Station  852 
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Figure  7-45.  Observed  and  computed  dissolved  oxygen  profiles  at  Station  852  for  1995 
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Figure  7-46.  Observed  and  computed  dissolved  oxygen  profiles  at  Station  852  for  1996 
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Figure  7-47.  Observed  and  computed  dissolved  oxygen  profiles  at  Station  852  for  1997 
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Figure  7-50.  Scatterplot  of  observed  and  computed 
chlorophyll  at  Lake  Washington  stations 
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Figure  7-51.  Scatterplot  of  observed  and  computed  total 
nitrogen  at  Lake  Washington  stations 


Chapter  7  Observation-based  Calibration 


Nitrate-Nitrite  (mg/L) 

Lake  Washington  Stations 


Figure  7-52.  Scatterplot  of  observed  and  computed 
nitrate  at  Lake  Washington  stations 


Figure  7-53.  Scatterplot  of  observed  and  computed 
ammonium  at  Lake  Washington  stations 
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Figure  7-54.  Scatterplot  of  observed  and  computed  total 
phosphorus  at  Lake  Washington  stations 


Figure  7-55.  Scatterplot  of  observed  and  computed  ortho¬ 
phosphorus  at  Lake  Washington  stations 
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7-56.  Scatterplot  of  observed  and  computed 

dissolved  oxygen  at  Lake  Washington  stations 
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Figure  7-57.  Relative  error  in  three  model  applications 
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Introduction 

A  rich  database  of  process  measurements  exists  for  Lake  Washington.  These 
include  measures  of  primary  production,  respiration,  deposition,  and  sediment- 
water  fluxes.  Extensive  nutrient  budgeting  for  the  lake  has  also  been  conducted. 
None  of  these  analyses  corresponds  with  the  model  application  period.  Still,  it  is 
worthwhile  to  examine  model  performance  in  light  of  these  measures  and 
analyses.  The  majority  of  measures  were  conducted  near  present  Station  852 
(Figure  8-1).  This  is  the  deepest  portion  of  the  lake  and  is  removed  from  major 
riverine  inflows.  Unless  otherwise  noted,  model  results  are  from  the  column  of 
cells  corresponding  to  Station  852. 


Primary  Production 

Net  primary  production  was  measured  during  the  mid-1 970’s  (Devol  and 
Packard  1978;  Richey  1979)  and  again  in  1980  (Quay  et  al.  1986).  For 
comparison  with  the  model,  monthly-average  values  were  computed  from 
individual  measures.  Averages  comprised  as  few  as  2  (December)  to  as  many  as 
13  observations  (May).  Observations  indicate  two  peaks  in  production  (Figure  8- 
2).  The  first  occurs  in  April,  concurrent  with  the  spring  phytoplankton  bloom. 

The  second,  higher  peak  occurs  in  July  when  observed  chlorophyll  is  much  lower 
than  the  spring  maximum.  This  pattern  resembles  that  for  estuaries,  e.g., 
Chesapeake  Bay,  in  which  maximum  biomass  and  maximum  production  are  out 
of  phase  (Harding  et  al.  2002).  The  model  indicates  a  system  in  which  biomass 
and  production  are  in  phase.  Peak  production  occurs  concurrent  with  the  spring 
bloom.  Computed  production  corresponds  well  with  observed  mean  values 
during  the  spring  months  but  is  closer  to  the  observed  minimum  values 
otherwise.  The  mean  annual  observed  net  production  is  ~0.46  g  C  m 2  day  1  while 
the  model  production  is  0.25  g  C  m'2  day'1. 

Considerable  effort  has  been  devoted  to  successfully  modeling  primary 
production  in  estuarine  systems  (Cerco  and  Noel  in  preparation).  The  estuarine 
model  has  two  key  differences  from  the  Lake  Washington  application.  The  first 
is  the  use  of  a  quadratic  rather  than  linear  phytoplankton  predation  term.  The 
second  is  the  distinction  of  spring  and  summer  algal  groups.  An  additional 
difference  in  the  two  efforts  is  the  availability  of  an  extensive,  contemporary 
database  of  primary  production  observations  for  Chesapeake  Bay .  With  some 
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effort,  the  observed  Lake  Washington  production  pattern  can  be  modeled. 
Complete  confidence  in  modeled  production  requires,  however,  a  comprehensive 
consistent  database  of  contemporary  observations  for  comparison  with  the  model 


Figure  8-1.  Lake  Washington  sample  stations 
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Years 


Figure  8-2.  Observed  and  computed  net  primary  production.  Observations  are 
monthly  means  of  reported  values.  Daily  modeled  values  are  from 
Station  852 


Respiration 

Respiration  was  measured  on  six  dates,  at  multiple  depths,  in  1980  (Quay  et 
al.  1986).  Observations  near  the  surface  (Figure  8-3)  indicate  highest  respiration 
rates  occur  in  June  and  July.  Near  the  bottom,  little  seasonality  is  apparent 
(Figure  8-4),  and  average  respiration  is  much  less  than  at  the  surface  (Figure  8- 
5).  The  model  demonstrates  remarkable  correspondence  with  the  majority  of 
surface  respiration  measures.  Near  the  bottom,  computed  respiration  exceeds 
observed.  Both  model  and  observations  indicate  a  surface-to-bottom  gradient  in 
respiration  but  the  observations  fall  off  much  more  rapidly,  below  10  m,  than  the 
model. 
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Figure  8-5.  Mean  observed  and  computed  respiration  at  multiple  depths 


The  observations  indicate  lower  microbial  activity  at  great  depths  than  near 
the  surface.  The  lower  activity  may  be  due  partly  to  lower  temperature  near  the 
bottom.  Attempts  to  reproduce  the  observations  by  manipulating  temperature 
coefficients  in  the  model  were  not  successful,  however.  It  is  suggested  that  the 
dissolved  organic  carbon  found  below  the  thermocline  is  more  refractory  than  the 
dissolved  organic  carbon  in  surface  waters.  Two  approaches  can  be  adopted  to 
improve  representation  of  the  vertical  variation  in  respiration.  The  first  is  to 
differentiate  labile  and  refractory  dissolved  organic  carbon  in  the  model.  This 
approach  requires  substantial  activity  in  coding  and  also  requires  a  means  to 
distinguish  loads  and  boundary  conditions  for  the  multiple  carbon  forms.  A 
second,  more  empirical,  approach  is  to  employ  lower  first-order  dissolved 
organic  carbon  mineralization  rates  below  the  thermocline. 


Deposition  and  Sediment-Water  Fluxes 

A  set  of  observed  deposition  and  sediment-water  fluxes  was  assembled  from 
multiple  sources.  These  were  converted  to  model  units  and  reported  in  Table  8-1 . 
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Table  8-1 

Observed  Deposition  and  Sediment-Water  Fluxes 

Quantity 

Range 

Comment 

Source 

Sediment  oxygen  demand, 
g  m*2  day'1 

0.05  to 

0.11 

Derived  from  diffusive 

C02  flux 

Kuivila  and  Murray 
(1984) 

Sediment  oxygen  demand, 
g  m'2  day'1 

0.14  to 

0.19 

Measured  oxygen 
consumption 

Kuivila  et  al.  (1988) 

Sediment  NH4  release, 
mg  N  m'2  day'1 

I. 5  to 

II. 4 

Diffusive  flux 

Kuivila  and  Murray 
(1 984) 

Sediment  P04  release, 
mg  P  m'2  day'1 

1.3 

Derived  from  P  budget 

Edmondson  and 

Lehman  (1981) 

Carbon  deposition, 
g  C  m'2  day'1 

0.15 

Average  rate  near 
bottom 

Kuivila  and  Murray 
(1984) 

Carbon  deposition, 
g  C  m’2  day'1 

0.05  to 

0.3 

Sediment  trap  at  50  m 

Quay  et  al.  (1986) 

Nitrogen  deposition, 
mg  N  m'2  day’1 

8  to  33 

Sediment  traps  at 
various  depths 

Edmondson  and 

Lehman  (1981) 

Phosphorus  deposition, 
mg  P  m*2  day'1 

2  to  67' 

Sediment  traps  at 
various  depths 

Edmondson  and 

Lehman  (1981) 

The  model  demonstrates  a  system  in  which  maximum  carbon  deposition 
occurs  in  June  and  July  (Figure  8-6).  The  deposition  peak  is  sharp  and  lags  by 
several  months  the  peak  in  net  production.  The  model  behavior  suggests  the 
maximum  deposition  is  a  result  of  the  high  first-order  phytoplankton  predation 
term  employed  in  summer.  Model  deposition  represents,  in  effect,  deposition  of 
zooplankton  fecal  matter.  The  peak  modeled  deposition  exceeds  deposition  rates 
measured  in  particle  traps,  but  the  model  mean,  0.26  g  C  m'2  day'1,  is 
characteristic  of  observed  deposition  rates.  Model  sediment  oxygen  demand 
occupies  a  range  from  0.2  to  0.4  g  m'2  day’1 .  Modeled  demand  is  at  the  upper 
range  of  observed  demand.  Most  of  the  observations  were  collected  in  colder 
months  (February-April  and  October)  during  which  temperature  effects  would  be 
expected  to  diminish  sediment  oxygen  demand.  In  view  of  the  sporadic  nature  of 
the  observations,  it  is  believed  that  the  differences  between  computed  and 
observed  sediment  oxygen  demand  are  insignificant. 

Modeled  nitrogen  deposition  follows  the  pattern  of  carbon  deposition  (Figure 
8-7).  As  with  carbon,  peak  modeled  deposition  exceeds  observed,  but  the  mean 
modeled  rate,  3 1  mg  N  m'2  day’1,  is  representative  of  the  observations.  Modeled 
sediment  ammonium  release  reflects  very1  well  the  observed  values.  A  novel 
result  from  the  model  is  the  indication  that  sediment  nitrate  uptake  is  a  substantial 
component  of  the  sediment-water  dissolved  nitrogen  flux.  No  observations  exist 
to  verify  the  model  results,  but  the  computed  values  are  expected  and  reasonable 
based  on  observations  conducted  in  estuaries  and  tidal  fresh  water  (DiToro 
2001). 


8-6 


Chapter  8  Process-based  Verification 


Figure  8-6.  Computed  carbon  deposition  and  sediment  oxygen  demand  at  Station 
852 


Figure  8-7.  Computed  particulate  and  dissolved  sediment-water  nitrogen  flux  at 
Station  852 


In  contrast  to  carbon  and  nitrogen,  peak  modeled  phosphorus  deposition 
(Figure  8-8)  is  less  than  peak  observed.  The  modeled  deposition  is  well  within 
the  observed  range.  Computed  phosphorus  flux  is  usually  representative  of  the 
observed  flux.  Larger  spikes  are  occasionally  computed  by  the  model.  These 
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occur  when  computed  dissolved  oxygen  concentration  dips  below  the  value  used 
in  the  sediment  model,  2  g  m'3,  to  initiate  transformation  to  anoxic  conditions. 
The  model  is  believed  to  successfully  represent  phosphorus  deposition  and 
sediment  release  in  Lake  Washington. 


Figure  8-8.  Computed  particulate  and  dissolved  sediment-water  phosphorus  flux 
at  Station  852 


System  Nutrient  Budgets 

Systemwide  budgets  were  constructed  for  nitrogen  (Figure  8-9)  and 
phosphorus  (Figure  8-10)  in  Lake  Washington  (Lake  Union  was  not  considered). 
Computed  sources  and  sinks  for  three  years  were  averaged  into  monthly  values. 

The  budgets  indicate  a  system  in  which  largest  nutrient  loads  occur  from 
November  through  March.  These  loads  are  balanced  by  outflow  through  Lake 
Union  and  by  net  accumulation  in  the  water  column.  Peak  nitrogen  deposition  to 
the  sediments  occurs  from  May  through  July.  The  nitrogen  deposition  is  balanced 
by  a  decrease  in  water  column  mass.  Dissolved  nitrogen  flux  is  into  the 
sediments  and  is  a  small  part  of  the  total  budget.  The  net  dissolved  nitrogen  flux 
into  the  sediments  is  the  result  of  sediment  nitrate  uptake  that  exceeds  sediment 
ammonium  release.  Peak  phosphorus  deposition  to  the  sediments  occurs  from 
April  to  July.  This  deposition  is  balanced  by  recycling  of  dissolved  phosphorus 
from  sediments  to  water  and  by  a  decrease  in  water  column  mass. 
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8-9.  Lake  Washington  model  nitrogen  budget.  Positive  fluxes  are  to  the 
Lake  Washington  water  column,  pnlws  =  particulate  nitrogen  flux  from 
water  to  sediment;  dnfws  =  dissolved  nitrogen  flux  at  sediment-water 
interface;  tn_flux  =  transport  at  Lake  Union  interface;  load  =  sum  of  all 
external  loads;  and  net  accumulation  =  change  in  water  column  mass 


Figure  8-10.  Lake  Washington  model  phosphorus  budget.  Positive  fluxes  are  to 
the  Lake  Washington  water  column,  ppfws  =  particulate  phosphorus 
flux  from  water  to  sediment;  benpo4  =  dissolved  phosphorus  flux  at 
sediment-water  interface;  tp_flux  =  transport  at  Lake  Union  interface; 
load  =  sum  of  all  external  loads;  and  net  accumulation  =  change  in 
water  column  mass 
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On  an  annual  basis,  the  Sammamish  River  is  computed  as  the  largest  source 
of  nutrients  to  Lake  Washington,  followed  by  the  Cedar  River  and  other 
distributed  sources  (Table  8-2).  The  majority  of  the  nutrient  load  is  deposited  in 
the  sediments.  A  lesser  amount  leaves  via  Lake  Union.  The  authors  compared 
their  budget  to  the  average  of  three  most  recent  years,  1976-1978,  reported  by 
Edmondson  and  Lehman  (1981)  (Table  8-3).  The  nutrient  loads  of  the  present 
study  are  30  percent  (nitrogen)  to  60  percent  (phosphorus)  higher  than  the  loads 
from  the  late  1970’s.  To  some  extent,  the  higher  loads  can  be  attributed  to  higher 
flows.  Mean  hydrologic  loading  reported  by  Edmondson  and  Lehman  for  1976- 
1978  was  33  m3  s'1  versus  39  m3  s'1  for  this  study.  The  18  percent  higher  flows 
during  this  study  do  not  completely  account  for  the  higher  loads.  It  cannot  be  said 
with  certainty  whether  the  remaining  disparity  indicates  a  loading  trend  in  the 
lake  or  differences  in  computational  methods.  Both  this  model  and  Edmondson 
and  Lehman  indicate  half  the  fluvial  nitrogen  load  is  deposited  in  the  sediments. 
Since  the  model  load  is  larger,  the  model  deposition  exceeds  the  Edmondson  and 
Lehman  budget.  Both  budgets  also  show  close  correspondence  in  the  fraction  of 
fluvial  phosphorus  load,  ~80  percent  deposited  in  the  sediments.  As  with 
nitrogen,  the  larger  phosphorus  load  of  this  study  produces  a  larger  net  load  to 
the  sediments.  Overall,  it  is  concluded  that  the  nutrient  budget  of  this  study 
provides  excellent  correspondence  with  the  independent  calculations  of 
Edmondson  and  Lehman. 


Table  8-2 

Model  Nutrient  Budget,  kg 

day’1 

Source 

Nitrogen 

Phosphorus 

Cedar 

861 

63.5 

Sammamish 

1251 

75.9 

Thornton 

55 

3.4 

Distributed 

756 

38.3 

CSO 

16 

2.3 

Atmospheric 

95 

24.7 

Particle  Flux 

-1271 

-198.2 

Dissolved  Flux 

-300 

53.8 

Outflow 

-1177 

-55.6 

Table  8-3 

Model  Nutrient  Budget  Compared  to  Edmondson  and  Lehman 
(1981) 

Source 

Nitrogen,  kg  day'1 

Phosphorus,  kg  day'1  j 

Model 

E&L 

Model 

E&L 

Fluvial 

2938 

2249 

183 

115 

Atmospheric 

95 

89 

25 

23 

Outflow 

1177 

1149 

56 

47 

Net  Input 

1571 

1190 

144 

92 
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9  Summary,  Conclusions, 
Recommendations 


Introduction 

In  September  2000,  the  King  County  Department  of  Natural  Resources  and 
Parks,  Washington  State,  signed  a  cooperative  agreement  with  the  U.S.  Army 
Engineer  Research  and  Development  Center  (ERDC)  to  provide  the  water  quality 
models  necessary  for  the  Sammamish-Washington  Analysis  Program.  ERDC 
responsibilities  under  the  agreement  included  the  following: 

1 .  Configure  and  execute  a  model  for  Lake  Washington. 

2.  Submit  written  project  reports. 

3.  Provide  training  and  expertise  to  county  personnel. 

4.  Assist  in  assessing  model  operation,  interpretation  of  field  data,  and 
model  results. 

5.  Modify  code  as  required. 

6.  Provide  expertise  on  lake  and  river  dynamics  such  as  nutrient  dynamics, 
ecology,  and  sediment-water  interactions. 

For  the  modeling  effort,  ERDC  applied  the  CH3D-WES  hydrodynamic 
model  and  the  CE-QUAL-ICM  eutrophication  model.  The  models  were  applied 
to  Lakes  Washington  and  Union  for  the  period  1995-1997.  Emphasis  was  on 
Lake  Washington.  Lake  Union  was  modeled  so  that  the  downstream  boundary 
condition  could  be  located  away  from  the  region  of  greatest  interest.  This  report 
presents  the  primary  documentation  of  the  model  effort. 


Flows,  Loads,  Boundary  Conditions 

Flows  into  Lake  Washington  may  be  divided  into  two  categories:  gauged  and 
distributed.  Gauged  flows  include  the  Sammamish  River,  the  Cedar  River,  and 
Thornton  Creek.  The  distributed  flows  enter  at  undefined  locations  and  are  not 
monitored.  These  were  computed  by  the  sponsor  based  on  precipitation  and  basin 
characteristics.  Daily,  gauged  flows  were  input  to  the  model  at  corresponding 
locations.  The  distributed  flows  were  pooled  into  10  locations  that  could  be 
identified  as  named  tributaries  or  features  and  input  to  the  model  on  a  daily  basis. 
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For  the  study  period,  the  Cedar  River  provided  60  percent  of  the  runoff  followed 
by  the  Sammamish  River  with  28  percent.  Thornton  Creek  provided  «1  percent. 
Distributed  flows  contributed  the  remaining  1 1  percent. 

Daily  runoff  loads  were  computed  as  the  product  of  flow  and  concentration. 
Concentrations  in  gauged  flows,  observed  at  approximately  monthly  intervals, 
were  provided  by  the  sponsor.  Linear  interpolation  was  used  to  create  daily 
concentrations  from  the  monthly  values.  Combined-sewer  overflows  (CSOs), 
sanitary-sewer  overflows,  and  spills  were  reported  at  28  locations.  Daily  flow 
volume  for  each  incident  was  provided  by  the  sponsor.  The  sponsor  also 
provided  characteristic  concentrations  for  the  discharges.  Daily  loads  were 
computed  as  the  product  of  flow  and  concentration  and  input  to  the  model  at 
appropriate  locations.  Annual  atmospheric  loads  were  obtained  from  Edmondson 
and  Lehman  (1981)  and  converted  to  daily  loads  for  input  to  the  model. 

The  Sammamish  River  was  the  predominant  nutrient  source  even  though 
flow  was  higher  in  the  Cedar.  The  Sammamish  contributed  greater  loads  because 
this  river  was  nutrient-enriched  relative  to  the  Cedar.  Nutrient  concentrations  in 
distributed  loads,  as  exemplified  by  Thornton  Creek,  were  enriched  relative  to 
both  large  rivers.  Although  CSO  events  contributed  enormous  loads  over  short 
periods,  these  sources  were  minor  compared  with  runoff  on  an  annual  basis. 
Likewise,  atmospheric  deposition  was  a  minor  component  of  the  annual  nutrient 
budget. 

On  an  annual  basis,  major  fecal  coliform  loads  were  from  the  Sammamish 
and  Cedar,  followed  by  distributed  loads.  CSO  loads  included  several  events  that 
rivaled  loads  from  major  runoff  sources,  but  spills  and  overflows  were  minor 
sources  on  a  long-term  basis. 

Suspended  solids  concentrations  in  the  Cedar  exceeded  the  Sammamish  so 
that  the  Cedar  was  the  dominant  source  to  the  system.  Distributed  sources  were 
enriched  with  solids  relative  to  major  runoff  sources,  but  the  role  of  distributed 
solids  sources  was  less  than  distributed  nutrient  sources.  Sewer  overflows  and 
spills  were  a  minor  component  of  the  long-term  suspended  solids  budget. 


Water  Quality  Model  Formulation 

The  CE-QUAL-ICM  water  quality'  model  was  applied  to  Lake  Washington. 
The  foundation  of  CE-QUAL-ICM  is  the  solution  to  the  three-dimensional  mass- 
conservation  equation  for  a  control  volume.  Control  volumes  correspond  to  cells 
on  the  model  grid.  As  applied  to  Lake  Washington,  the  CE-QUAL-ICM  model 
incorporated  20  state  variables  in  the  water  column  including  physical  variables, 
phytoplankton,  multiple  forms  of  carbon,  nitrogen,  and  phosphorus,  a  pathogen, 
and  two  toxicants  (Table  9-1).  The  pathogen  and  toxicants  were  added 
specifically  for  the  Lake  Washington  application. 


9-2 


Chapter  9  Summary,  Conclusions,  Recommendations 


Table  9-1 

Lake  Washington  Water  Quality  Model  State  Variables 

Temperature 

Fixed  Solids 

Phytoplankton 

Fecal  Coliform 

Dissolved  Organic  Carbon 

Labile  Particulate  Organic  Carbon 

Refractory  Particulate  Organic  Carbon 

Ammonium 

Nitrate-*-  Nitrite 

Dissolved  Organic  Nitrogen 

Labile  Particulate  Organic  Nitrogen 

Refractory  Particulate  Organic  Nitrogen 

Total  Phosphate 

Dissolved  Organic  Phosphorus 

Labile  Particulate  Organic  Phosphorus 

Refractory  Particulate  Organic  Phosphorus 

Chemical  Oxygen  Demand 

Dissolved  Oxygen 

Toxicant  One 

Toxicant  Two 

Modeling  Processes  at  the  Sediment-Water 
Interface 

The  predictive  benthic  sediment  model  applied  to  Lake  Washington  was  first 
developed  for  use  in  Chesapeake  Bay  (DiToro  2001).  The  model  (Table  9-2)  is 
driven  by  net  settling  of  organic  matter  from  the  water  column  to  the  sediments. 
In  the  sediments,  the  model  simulates  the  diagenesis  (decay)  of  the  organic 
matter.  Diagenesis  produces  oxygen  demand  and  inorganic  nutrients.  Oxygen 
demand,  as  sulfide  (in  salt  water)  or  methane  (in  fresh  water),  takes  three  paths 
out  of  the  sediments:  export  to  the  water  column  as  chemical  oxygen  demand, 
oxidation  at  the  sediment-water  interface  as  sediment  oxygen  demand,  or  burial 
to  deep,  inactive  sediments.  Inorganic  nutrients  produced  by  diagenesis  take  two 
paths  out  of  the  sediments:  release  to  the  water  column  or  burial  to  deep,  inactive 
sediments. 


Table  9-2 

Sediment  Model  State  Variables  and  Fluxes 

State  Variable 

Sediment-Water  Flux 

Temperature 

Particulate  Organic  Carbon 

Sediment  Oxygen  Demand 

Sulfide/Methane 

Release  of  Chemical  Oxygen  Demand 

Particulate  Organic  Nitrogen 

Ammonium 

Ammonium  Flux 

Nitrate 

Nitrate  Flux 

Particulate  Organic  Phosphorus 

Phosphate 

Phosphate  Flux 

Coupling  with  the  Hydrodynamic  Model 

Transport  processes  were  modeled  by  a  three-dimensional  hydrodynamic 
model  that  operated  independently  of  the  water  quality  model  (Johnson  et  al. 


Chapter  9  Summary,  Conclusions,  Recommendations 


9-3 


2003).  Transport  infonnation  from  the  hydrodynamic  model  was  processed  and 
stored  on-line  for  subsequent  use  by  the  water  quality  model.  Computed  flows, 
surface  elevations,  and  vertical  diffusivities  were  output  at  1-hour  intervals  for 
use  in  the  water  quality  model.  The  algorithms  and  codes  for  linking  the 
hydrodynamic  and  water  quality  models  were  developed  over  a  decade  ago  and 
have  been  tested  and  proved  in  extensive  applications  since  then.  For  this 
application,  linkage  was  verified  through  comparison  of  temperatures  computed 
by  die  two  models.  A  total  of  3,468  comparisons  of  temperatures  computed  at  18 
stations  over  a  3-year  period  were  made.  Correspondence  between  the  two 
models  was  excellent.  The  linear  correlation  coefficient  of  the  independently 
computed  temperatures  was  -0.99  and  relative  error  was  less  than  1  percent. 

A  great  deal  of  difficulty  was  encountered  specifying  hydrodynamic 
boundary  conditions  at  the  lock  and  dam  on  the  lower  end  of  Lake  Union.  After 
much  trial  and  error,  outflow  was  specified  based  on  volume  conservation.  This 
specification  was  sufficient  for  the  Lake  Washington  application.  Modeling  of 
Lake  Union  will  require  much  more  attention  to  flow  and  salinity  exchange  at  the 
structure.  Both  the  structure  itself  and  the  operation  of  the  lock  will  have  to  be 
considered.  Modeling  of  salinity  intrusion  and,  consequently,  circulation  in  Lake 
Union  will  be  impossible  without  improvement  in  boundary  condition 
specification.  One  recent  study  (Johnson  et  al.  2004)  may  provide  some  insight 
into  the  problem.  The  challenge  of  modeling  this  boundary  should  not  be 
underestimated. 


Observation-based  Calibration 

The  primary  database  for  calibration  of  the  model  was  derived  from  a 
monitoring  program  conducted  by  King  County  Department  of  Parks  and  Natural 
Resources.  The  model  was  compared  to  observed  temperature,  dissolved  oxygen, 
ammonium,  nitrate+nitrite,  total  nitrogen,  ortho-phosphorus,  total  phosphorus, 
chlorophyll,  disk  visibility,  and  fecal  coliform.  Comparisons  were  produced  in  a 
variety  of  graphical  and  statistical  formats.  Graphical  formats  included 
examination  of  time  series,  vertical  profiles,  and  scatterplots.  Statistics  computed 
included  mean  error,  relative  error,  and  variance. 

Overall,  the  model  provided  excellent  representation  of  the  annual  cycle  of 
temperature,  chlorophyll,  nutrients,  and  dissolved  oxygen  in  the  lake.  The 
vertical  structure  was  also  well-represented.  In  relative  terms,  the  greatest 
disparities  between  computations  and  observations  were  for  ortho-phosphorus 
and  total  phosphorus.  Least  relative  error  was  for  ammonium,  followed  by 
chlorophyll  and  dissolved  oxygen,  which  had  nearly  equivalent  relative  error. 
Relative  error  in  the  Lake  Washington  application  compared  most  favorably  with 
several  alternate  model  applications. 

One  problem  area  in  the  model  application  was  the  trapping  of  excess  ortho¬ 
phosphorus  below  the  thermocline  from  late  spring  through  fall  turnover.  The 
computation  of  excess  ortho-phosphorus  has  been  experienced  in  alternate  model 
applications.  The  excess  is  believed  to  originate  in  the  lack  of  bacterial 
phosphorus  uptake  in  the  model.  It  is  not  known  that  any  models  presently 
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include  this  process.  Mechanistic  modeling  of  bacterial  activity  is  probably 
unwarranted,  but  some  empirical  representation  of  bacterial  uptake  is  possible 
and  should  be  attempted.  The  sponsor  has  suggested  that  binding  of  phosphorus 
in  zooplankton  fecal  pellets,  and  subsequent  rapid  sinking  of  these  pellets,  is  not 
represented  in  the  present  model.  Phosphorus  removal  via  sinking  may  remove 
phosphorus  that  is  presently  being  mineralized  and  causing  an  excess.  This 
possibility  is  worth  examining.  Since  ortho-phosphorus  is  well -represented  in 
surface  waters  of  the  lake,  where  utilization  as  an  algal  nutrient  is  important,  the 
computed  excess  of  ortho-phosphorus  does  not  compromise  the  utility  of  the 
model.  However,  the  problem  should  be  addressed. 

Fecal  coliform  bacteria  are  indicators  of  pollution  from  human  and/or  animal 
waste.  This  study  identifies  the  sources  of  coliform  pollution  and  computes 
concentrations  consistent  with  observations.  The  highest  observed  values  are  not 
represented  in  the  model.  The  computations  are  affected  by  the  crude  loading 
information;  e.g.,  all  CSO  events  are  characterized  by  a  single  concentration. 
Improved  calculations  of  fecal  coliform  will  require  improved  loading  estimates. 

The  measure  of  light  attenuation  available  for  this  study  was  disk  visibility, 
also  known  as  Secchi  depth.  The  corresponding  model  variable  is  diffuse  light 
attenuation.  Attenuation  can  be  related  to  disk  visibility  through  an  empirical 
relationship.  The  model  properly  represents  the  magnitude  of  disk  visibility  and 
the  processes  that  contribute  to  seasonal  variations  in  this  observed  quantity.  The 
observations  show  a  great  deal  of  scatter  about  the  model  value  and  a  great  deal 
of  scatter  in  general.  The  loose  agreement  between  model  and  observed  visibility 
is  attributed  to  three  factors:  the  inherent  uncertainty  in  the  disk  visibility 
measure,  lack  of  solids  data  for  calibration  of  this  attenuating  factor,  and  the 
conversion  from  attenuation  to  disk  visibility.  Consideration  should  be  given  to 
replacement  of  the  disk  visibility  measure  with  direct  measurement  of  subsurface 
irradiance,  from  which  attenuation  can  be  readily  calculated.  The  sponsor  should 
consider  measurement  of  color  and  turbidity  and  the  relationship  of  these  to 
attenuation.  The  present  model  uses  an  additive  formula  that  computes 
attenuation  as  a  function  of  background  attenuation,  suspended  solids,  and 
chlorophyll.  If  sufficient  observations  become  available,  consideration  should  be 
given  to  application  of  an  optical  model  based  on  absorption  and  scattering  (e.g., 
Gallegos  et  al.  1990). 


Process-based  Verification 

A  rich  database  of  process  measurements  exists  for  Lake  Washington.  These 
include  measures  of  primary  production,  respiration,  deposition,  and  sediment- 
water  fluxes.  Extensive  nutrient  budgeting  for  the  lake  has  also  been  conducted. 
Although  none  of  these  analyses  corresponds  with  the  model  application  period, 
examination  of  model  performance  in  light  of  these  measures  and  analyses  was 
thought  worthwhile. 

Observations  indicate  two  peaks  in  net  primary  production.  The  first  occurs 
in  April,  concurrent  with  the  spring  phytoplankton  bloom.  The  second,  higher 
peak  occurs  in  July  when  observed  chlorophyll  is  much  lower  than  the  spring 
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maximum.  This  pattern  resembles  estuaries,  e.g.,  Chesapeake  Bay,  in  which 
maximum  biomass  and  maximum  production  are  out  of  phase.  The  model 
indicates  a  system  in  which  biomass  and  production  are  in  phase.  Peak 
production  occurs  concurrent  with  the  spring  bloom.  Computed  production 
corresponds  well  with  observed  mean  values  during  the  spring  months  but  is 
closer  to  the  observed  minimum  values  otherwise.  The  mean  annual  observed  net 
production  is  -0.46  g  C  m'2  day'1  while  the  model  production  is  0.25  g  C  m'2 
day'1. 

Observations  near  the  surface  indicate  highest  respiration  rates  occur  in  June 
and  July.  Near  the  bottom,  little  seasonality  is  apparent  and  average  respiration  is 
much  less  than  at  the  surface.  The  model  demonstrates  remarkable 
correspondence  with  the  majority  of  surface  respiration  measures.  Both  model 
and  observations  indicate  a  surface-to-bottom  gradient  in  respiration,  but  the 
observations  fall  off  much  more  rapidly,  below  10  m,  than  the  model. 

Modeled  particulate  carbon,  nitrogen,  and  phosphorus  deposition  to  the 
bottom  sediments  are  characteristic  of  rates  observed  in  sediment  traps.  Modeled 
sediment  oxygen  demand  is  at  the  upper  range  of  observed  demand.  Most  of  the 
observations  were  collected  in  colder  months  during  which  temperature  effects 
diminish  sediment  oxygen  demand.  In  view  of  the  sporadic  nature  of  the 
observations,  the  differences  between  computed  and  observed  sediment  oxygen 
demand  are  believed  to  be  insignificant.  Modeled  sediment  ammonium  release 
reflects  very  well  the  observed  values.  A  novel  result  from  the  model  is  the 
indication  that  sediment  nitrate  uptake  is  a  substantial  component  of  the 
sediment-water  dissolved  nitrogen  flux.  Computed  phosphorus  flux  is  usually 
representative  of  the  observed  flux. 

On  an  annual  basis,  the  Sammamish  River  is  computed  to  be  the  largest 
source  of  nutrients  to  Lake  Washington,  followed  by  the  Cedar  River  and  other 
distributed  sources.  The  majority  of  the  nutrient  load  is  deposited  in  the 
sediments.  A  lesser  amount  leaves  via  Lake  Union.  The  budget  in  this  study  was 
compared  to  the  average  of  three  most  recent  years,  1976-1978,  reported  by 
Edmondson  and  Lehman  (1981).  The  nutrient  loads  of  this  study  are  30  percent 
(nitrogen)  to  60  percent  (phosphorus)  higher  than  the  loads  from  the  late  1970’s. 

It  cannot  be  said  with  certainty  whether  the  disparity  indicates  a  loading  trend  in 
the  lake  or  differences  in  computational  methods. 

Considerable  effort  has  been  devoted  to  successfully  modeling  primary 
production  in  estuarine  systems.  The  present  estuarine  model  (Cerco  and  Noel  in 
preparation)  has  two  algal  groups  and  two  zooplankton  groups.  The  estuarine 
model  uses  a  quadratic  term  to  represent  additional  predation  on  phytoplankton 
by  microheterotrophs.  The  quadratic  term  has  been  found  to  be  a  key  to 
computing  both  phytoplankton  biomass  and  production.  It  is  understood  that  the 
addition  of  zooplankton  to  the  Lake  Washington  model  is  planned.  The  addition 
of  multiple  algal  groups  is  also  under  consideration.  The  addition  of  zooplankton 
will  require  modification  of  the  present  algal  predation  term.  Adoption  of  a 
quadratic  formulation  is  recommended.  The  sponsor  may  also  wish  to  explore 
alternate  formulations.  With  some  effort,  the  observed  Lake  Washington 
production  pattern  can  be  modeled.  Complete  confidence  in  modeled  production 
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requires,  however,  a  comprehensive,  consistent  database  of  contemporary 
observations  for  comparison  with  the  model. 

Two  approaches  can  be  adopted  to  improve  representation  of  the  vertical 
variation  in  respiration.  The  first  is  to  differentiate  labile  and  refractory  dissolved 
organic  carbon  in  the  model.  This  approach  requires  substantial  activity  in  coding 
and  also  requires  a  means  to  distinguish  loads  and  boundary  conditions  for  the 
multiple  carbon  forms.  A  second  approach  is  to  employ  lower  first-order 
dissolved  organic  carbon  mineralization  rates  below  the  thermocline.  This  second 
approach  is  recommended. 

Lake  Washington  is  the  first  application  of  the  sediment  diagenesis  model  to 
a  lake.  It  is  believed  that  the  results  are  most  satisfactory.  The  successful 
application  is  due  partly  to  the  absence  of  anoxia  in  Lake  Washington.  Under 
anoxic  conditions,  the  diagenesis  model  represents  estuarine  conditions  in  which 
reduction  of  sulfate  to  sulfide  is  the  dominant  process.  Reduction  and  oxidation 
of  iron  are  not  addressed  in  the  present  model  code  although  developmental  work 
has  been  conducted  in  this  area  (DiToro  2001).  Application  of  the  model  to  an 
anoxic  lake,  such  as  Lake  Sammamish,  may  require  inclusion  of  iron  dynamics 
into  the  diagenesis  and  water  column  models. 

Apparent,  substantial  increases  in  nutrient  loading  to  Lake  Washington  since 
the  late  1970’s  were  detected.  Extensive  analyses  were  not  conducted  to  verify 
the  apparent  increase.  Certainly,  increased  nutrient  loading  is  consistent  with 
development  trends  in  the  watershed  and  with  increased  loading  observed  to 
coastal  systems  nationwide.  The  sponsor  should  conduct  necessary  sampling  and 
analyses  to  determine  if  nutrient  loading  to  Lake  Washington  is,  indeed, 
increasing. 
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The  CE-QUAL-ICM  three-dimensional  eutrophication  model  was  applied  to  Lake  Washington,  northwestern  Washington  State,  for  the 
period  1995-1997.  Transport  processes  were  obtained  from  the  companion  CH3D-WES  hydrodynamic  model.  The  model  activated  18 
state  variables  in  the  water  column  including  physical  variables;  phytoplankton;  multiple  forms  of  carbon,  nitrogen,  and  phosphorus, 
dissolved  oxygen;  and  fecal  cohform.  The  model  provided  excellent  representation  of  the  annual  cycle  of  temperature,  chlorophyll, 
nutrients,  and  dissolved  oxygen  in  the  lake.  The  vertical  structure  was  also  well-represented.  The  water  column  was  coupled  to  a 
predictive  sediment  diagenesis  model  that  computed  sediment -water  fluxes  of  dissolved  oxygen,  methane,  ammonium,  nitrate,  and 
phosphate,  based  on  computed  inputs  of  particulate  organic  matter.  Computed  values  of  deposition  and  fluxes  were  in  reasonable 
agreement  with  published  values.  Nutrient  loads  were  calculated  and  nutrient  budgets  were  constructed  as  part  of  the  model  exercise. 
Load  sources  included  river  inflows,  distributed  loads,  sewer  overflows,  and  atmospheric  loading.  The  Sammamish  River  was  identified 
as  the  largest  source  of  nutrients  to  Lake  Washington,  followed  by  the  Cedar  River  and  other  distributed  sources.  The  majority  of  the 
nutrient  load  is  deposited  in  the  sediments.  A  lesser  amount  leaves  via  Lake  Union.  Nutrient  loads  in  this  study  were  30  percent 
(nitrogen)  to  60  percent  (phosphorus)  higher  than  the  loads  from  the  late  1 970’s. 
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